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Summary 


This report describes the development of adaptive importance sampling tech¬ 
niques for estimating false alarm probabilities of detectors that use space-time 
adaptive processing (STAP) algorithms. The first paper to appear on this topic 
is [9]; it lays the groundwork for developing powerful estimation algorithms, 
based on the work in fast simulation carried out by the author of [5]. The re¬ 
port also presents some new detection algorithms that are intended to be robust 
under various conditions. 

Fast simulation using importance sampling methods have been notably suc¬ 
cessful in the study of conventional constant false alarm rate radar detectors, 
and in several other applications. The principal objectives here are to examine 
the viability of using importance sampling methods for STAP detection, develop 
these methods into powerful analysis and design algorithms, and use them for 
synthesizing novel detection structures. 

Several STAP detectors have been investigated from the standpoint of ap¬ 
plying importance sampling to characterize their performances. Various bias¬ 
ing techniques have also been devised and implemented, resulting in significant 
speed-ups in performance evaluation compared to conventional Monte Carlo 
methods. The important problem of detector threshold determination has been 
addressed and solved by fast simulation. 

Robust variants such as the envelope-law and geometric-mean detectors for 
STAP processing have been suggested, their CFAR property established, and 
performance thoroughly evaluated using IS techniques. It has been shown that 
their detection performances are decidedly better than those of their conven¬ 
tional square-law counterparts when training data are contaminated by interfer¬ 
es, while maintaining almost equal detection performances under homogeneous 
conditions. The work reported here paves the way to development of more 
advanced estimation techniques that can facilitate design of powerful and ro¬ 
bust detection algorithms designed to counter hostile and heterogenous clutter 
environments. 

During the pendency of the contract, which started in late 2004, the PI and 
his associates have presented their results at various conferences, in the UK, 
USA, Europe, and China. One journal paper 1 is due to appear in the IEEE 

1 A pre-print manuscript copy is attached at the end of this report. Conference presentations 
have not been included here. 



Transactions on AES in January 2007 and two more journal manuscripts are 
to be shortly submitted for review. Owing to unavoidable circumstances the 
contract has had to be terminated prematurely. Given the satisfactory results 
obtained in this short time frame, it is hoped that more work will be carried 
out along the same lines in the present research area. 
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Chapter 1 


Importance sampling for 
STAP detection analysis 
and design - the AMF 

1.1 Introduction 

Estimation of false alarm probabilities of detection algorithms that employ 
space-time processing is examined here using forced Monte Carlo or importance 
sampling simulation (IS). Space-time adaptive processing (STAP) algorithms 
are of much importance for radar detection. They are notoriously intensive 
from a computational point of view, with the more advanced (and robust) ones 
being also analytically difficult to quantize, [3]. Therefore it is appropriate to 
attempt to develop fast simulation methods that could be used in their analysis 
and design. 

In this chapter we use lessons learnt from developing IS techniques for char¬ 
acterizing conventional constant false alarm rate (CFAR) detectors, [5], and 
describe an experiment in applying them to STAP detection. The starting (and 
ending) point of this unpretentious effort is the celebrated adaptive matched 
filter (AMF) derived in [1] and which represents the arrayed version of the 
workhorse cell averaging CFAR detector for conventional radar signal process¬ 
ing algorithms. The false alarm probability (FAP) performance of the AMF 
detector is known in integral form and can be numerically computed to any 
desired accuracy. Thus it forms a suitable basis for validating our simulation 
experiments. Two specific IS methods (described in the sequel) are presented 
and the better (and also easier) one is implemented. On a general note, IS is 
the chief simulation methodology for rare event estimation. It is an enduring 
method that has distinguished itself in several areas of science and engineering. 
Briefly, IS works by biasing original probability distributions in ways that ac¬ 
celerate the occurrences of rare events, conducting simulations with these new 
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distributions, and then compensating the obtained results for the changes made. 
The principal consequence of this procedure is that unbiased probability esti¬ 
mates with low variances are obtained quickly. The main task in IS therefore 
is determination of good simulation distributions for an application, either as a 
one-shot feat or adaptively. Simulations performed using such distributions can 
provide enormous speed-ups if they are chosen with due care and mathematical 
precision. Indeed, if applied successfully, simulation lengths needed to estimate 
very low probabilities become (only) weakly dependent on the actual proba¬ 
bilities. It is thus possible to evaluate any probability in reasonable amounts 
of simulation time. There have been more recent attempts in the literature, 
for example [13], [14], to apply IS for FAP estimation of conventional CFAR 
detectors with varying degrees of success. 

During the conduct of simulations reported herein, some issues concerning 
the adaptive IS algorithms used arise, and these are discussed briefly. More in¬ 
vestigation is required into them. However, the positive outcome of the methods 
used is that excellent match with numerical results is obtained. The succeeding 
sections provide a short statement of the AMF algorithm, how IS biasing can 
be performed to hasten false alarm events, description of the so called ^-method 
which is a conditional IS technique developed originally for studying sums of 
random variables ([6]), the fast algorithms used, how inverse IS can be used 
to estimate (and therefore design) detector thresholds, simulation results, and 
a concluding discussion. The work carried out in this chapter has appeared in 
[9], 

1.2 The AMF detector 

In a radar system consisting of a linear array of N s antenna elements, a burst of 
N t pulses is transmitted and each element receives as many return samples in 
any one range gate. The N s N t = N samples are complex (because of I and Q 
channel processing) and are referred to as the primary data. They may contain 
a target and represent the range gate to be tested. The samples are arranged 
in an N x 1 column vector and denoted as x. The target return is modelled as 
consisting of a known direction vector s with an unknown complex amplitude in 
addition to clutter, interference, and noise. There are L other iV-length complex 
vectors, called the training data, obtained from as many nearby range gates and 

assumed to be free of target signal. These are denoted as x(Z), l = 1 ,_ ,L. 

It is assumed from now on that the training data is free of other targets or 
contamination 1 . The primary and secondary data vectors are assumed to be 
jointly independent and complex Gaussian, sharing the NxN covariance matrix 
R = £{XXt}, where the superscript f denotes complex conjugate transpose. 

x In the second chapter, detection performance is carried out for training data contaminated 
by interfering targets. 
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Under these assumptions the AMF detection test, as obtained in [1], is given by 


IstR-ixl 2 h , 

—- < V 

stR -1 S Ho 


( 1 . 1 ) 


where 

R= ^5Z J= i X (0 x (0 f 

is the estimated covariance matrix of x based on the training data (also referred 
to as sample matrix), and 77 is a threshold used to set the FAP at some desired 
level. This test has the CFAR property. The FAP a of the test is known to be 
given by 

L\ l' 1 x L ~ N+1 (l—x) N ~ 2 

(L — N + l)!(iV — 2 )! J Q (l + V x/L) L ~ N + 1 dX ( j 

which can be used to numerically determine the threshold setting for a desired 
FAP. As shown in [1], the test in (1.1) can be rewritten as 


IstR^xl 2 § 

7/S^R 1 S 

Ho 

= 

77S^R _1 RR~ 1 s 

= 

»?s t R _1 ^^ ;= i x (0 x (0 t R- _ 

= 

j- J2 l=1 s t R- _1 x(0x(0 t R _1 s 

= 



This is in the form of a vector (or, array) version of the usual CA-CFAR test. 
The LHS is a square law detector, being the output of a matched filter (matched 
to the direction s in which the array is steered) for incoherent detection using 
the so-called sample matrix inversion (SMI) beamformer weights R^ 1 s. The 
RHS represents a cell averaging term. Further details on these issues can be 
found in the references mentioned above. 


1.3 False alarm probability estimation using IS 

Two methods to quickly estimate FAPs are two-dimensional (2-d) biasing and 
the conditional (/-method procedure, described in this section. 

1.3.1 2-d biasing 

To estimate FAP using IS, we make the following observations. Suppose each 
complex sample of a training vector is scaled by a real number 6 1/<2 . This has 
the effect of scaling the covariance matrix estimate R by 9. Therefore, as far 
as the covariance estimate is concerned, both sides of the test in (1.3) remain 
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unaffected by the scaling. However, each training vector being scaled by 0 1 / 2 
results in a scaling of the RHS by 9. Hence choosing 9 less than unity will have 
the effect of compressing the density function of the random threshold of the 
test. Further, a scaling of each complex component of the primary vector by 
a real a 1 / 2 will achieve a scaling of the LHS of the test by a. Thus, choosing 
a larger and 9 smaller than unity will achieve an increase in the frequency 
of occurrence of a false alarm event during simulation. The IS optimization 
problem will be a two-parameter one. 

The (unbiased) IS estimator, using (1.1), can be expressed as 

S = lEfKIs^Xl 2 > ^s t R- 1 s)W(X,X L ;a,6l); ~ /* (1.4) 

where the notation ~ /* means that all random variables are drawn from biased 
distributions, and X L = (X(l),... ,X(L)) T with K denoting length of the IS 
simulation. In setting up their joint densities, we use the fact that the FAP 
of the AMF has the CFAR property and is independent of the true covariance 
matrix R. This is true under the assumption of Gaussian distributions for the 
data. In such a case, the simulation of the AMF test can be carried out for data 
possessing a diagonal covariance matrix I, denoting the N x N identity matrix. 
Therefore, primary and training data can be generated as complex vectors with 
independent components. The unbiased joint densities are 

e -x f x e -£f x(Z) f x(0 

/( x ) = —JT and /( x i) = -- Zln - 

7T 7T^ V 


so that 


/(x,x L )= 


g-x f x-£i x(J) + x(J) 


n (L+l)N 

With scaling performed as described above, the biased joint density takes the 
form 




0 — —X 1 X 


/*(X,X L ) n (L+l)N a NQLN 

resulting in the weighting function 

/( x , x l) 


W(X,X i; a,0) ^ 


/*( x ,x L ) 

= Ca N 9 LN e A / a e B ' e 


where 


A = x^x, B = x(Z)^x(Z), and C 

I l i variance of the IS estimator a can be expressed as 

vara = -t[/(iz) — a 2 ] 

A 


= P ~{A+B) 


(1.5) 


( 1 . 6 ) 
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where v is the vector biasing parameter ( a,9) T £ [l,oo) x (0,1], Denoting by 
A the false alarm event in (1.4), the quantity / is given by 

I(v) = E*{l{A)W 2 {X,X L -v)} 

= E{l(A)W(X,X L ;v)} (1.7) 

where the expectation If* proceeds over biased distributions. Minimization 
of var a with respect to the biasing parameters is equivalent to minimization 
of I and is described in the Appendix. Although not implemented here, this 
description has been included since it is foreseen that such a method could be 
useful in situations wherein the g-method might be difficult to apply. 

1.3.2 The g-method estimator 

This method exploits knowledge of underlying distributions more effectively, 
yielding a more powerful estimator. Additional advantages are that only a 
scalar parameter optimization problem needs to be tackled and the inverse IS 
problem (for threshold optimization or selection) can be easily solved. The FAP 
can be written as 

a = P(|s t R -1 X| 2 > r]s^B,~ 1 s\H 0 ) 

= ^{PflstR^XI 2 > rys t R- 1 s|X L ,R 0 )} 

= E{g(X L )} (1.8) 

Note that the conditioning in the second step above is equivalent to the condition 
that a covariance matrix estimate is given. We proceed to estimate a using the 
form in the third step above. 

With the condition in mind it is easy to show, assuming that X is rotationally 
invariant and Gaussian, that the random variable sTR _1 X = w^X is distributed 
as CAfi(0, w*Rw) with independent real and imaginary components, and the 
weight vector w = R _1 s. The random variable Y = |sfR -1 X| 2 therefore is 
exponential and has density function 

g-y/w 1 Rw 

/(y|Xi, Ho) = -——-, y > 0 

w'Kw 

Therefore 

g(X L ) = P(y ^jjstR-^IXi.Ro) 

g~ r] s*R" 1 s/w , Rw 

Note that if R = R, then g{Xi) = e~ v and this is the FAP of the AMF when 
the covariance matrix is known. As discussed before, we are simulating with 
homogeneous data possessing an identity covariance matrix, that is, with R = I. 
The g-method IS estimator then takes the form 

9{^l)W{X l -6) 
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e~ vD W{K L ', 0); ~ 


(1.9) 


1 

K 



where 


D 


stR-is 

|w| 2 

stR-is 

st(R-i)2 s 


( 1 . 10 ) 


Choosing the (single) biasing parameter 9 < 1 thus produces a decrease in D, 
thereby causing a higher frequency of occurrence of the false alarm event or, 
more appropriately in this case, a larger value of the (/-function. Note that use 
of the (/-method obviates the need to bias primary data vectors. Determination 
of a good value of 9 proceeds as before. The weighting function is simply 

W(x L ; 9) = e LN e- {1 - l,6)B (1.11) 


which can be deduced from (1.5) by setting a = 1. The scaling 9 is optimized 
by 

9 m+1 =9 m -5 e J^- ( 1 . 12 ) 

r g '(0m) 

which is just a one-dimensional version of (1.16). Estimates of the /-function 
and its derivatives are given by 

Ig{6) = ^E|V( x l)W 2 (X l ;0); 

I' g {0) = ^^V( x lW X l;W(X l ;0); ~/* 

%'(*) = 9 2 ^L)W(X L ;9)W ee (X L -,9y, ~ /* 

See appendix 1.5 for definition of the above quantities. 


Simulation gain 

A useful measure of the effectiveness of any IS procedure is the simulation 
gain T. It is the ratio of simulation lengths required by conventional MC and 
IS estimators to achieve the same error variance. Setting the variance in (1.6) 
equal to (a—a 2 )/k (being the MC variance) where k denotes the length required 
by the MC estimator, yields the gain 


I[y) - a 2 

While the simulation gain is useful in learning how much faster than MC an IS 
technique is in terms of simulation length, it also serves the purpose of comparing 
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different IS estimators. In actual simulations, an estimate of T is made by using 
the estimates for a and I. The (/-method estimator has simulation gain given 

by 

= a-a 2 
9 I g (0) -a 2 

where I g = E+{g 2 (XL)W 2 (XL',0)}, and it can be estimated during simula¬ 
tion. It always has a smaller variance and consequently larger gain than the 
IS estimator discussed in the previous section. Indeed, without IS (W = 1), 
I g = E{g 2 (l^i,)} < E{g(X.L)} = a. That I g < I with IS was proved in [5] for 
conventional CFAR detectors. The proof in the case of the detectors considered 
here is similar, and will be omitted. 

1.3.3 Inverse IS and threshold determination 

The inverse problem, namely that of finding by fast simulation the value of 
detector threshold 77 satisfying a prescribed FAP, is readily solved using the g- 
method estimator. This is done by minimizing the stochastic objective function 

J(v) = [a g (v) ~ Qo ] 2 

where a 0 is a desired FAP. An example is shown in Figure 1.1. It is clear that 
all detection algorithms that involve a threshold crossing will possess objective 
functions that have the general behaviour shown, assuming of course that the 
FAP estimate is a decreasing function of its argument g g . Minimization of J 
with respect to ?/ is carried out by the algorithm 

Vm+l — Vm H - Or/ — ^ •> — 1,2,... (1.13) 

a' g {rim) 

where S v is a step-size parameter and the derivative estimator in the denomi¬ 
nator is given by 


= De-^mx^ey, -/* (1.14) 

with the prime indicating derivative with respect to g. Note in passing that 
this derivative estimator actually estimates (negative of) the probability density 
function of the AMF statistic on the left hand side of (1.1) under H 0 . 

1.4 Numerical results 

The FAP a obtained by direct numerical integration of (1.2) is shown in Fig¬ 
ure 1.2 and is used for comparing IS results, which are displayed in the remaining 
figures. The AMF detector consists of L = 704 trainng vectors each of length 
N = 352. Shown in Figure 1.3 is one instance of adaptive IS estimation of FAP 
for a (known) threshold of rj = 56.50432. Figures 1.4 and 1.5 depict results from 
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Figure 1.1: An objective fucntion. 



Figure 1.2: Numerically computed FAP of the AMF detector. 


















Figure 1.3: Convergence of FAP using adaptive IS algorithms. 


AMF: L = 704, N = 352 



Figure 1.4: Threshold optimization for AMF detector using inverse IS. 
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Figure 1.5: FAP estimates resulting from threshold optimization algorithms. 


implementing the inverse IS algorithms. These are estimated threshold settings 
and associated FAPs respectively. It is evident that match with the results in 
Figure 1.2 is excellent and this has been numerically confirmed. 

Optimum biasing parameters are shown in Figure 1.6 and simulation gains 
obtained in Figure 1.7. Shown in Figure 1.8 is the gain as function of FAP. The 
number of trials K needed to provide an accuracy of ±10% with 95% confidence 
is shown in Figure 1.9. 

1.4.1 Discussion 

The IS simulation results obtained here appear deceptively smooth and certainly 
beg an obvious question. Indeed, an artifice has been employed here to generate 
them. It was used by the first author in previous work on capacity estimation of 
MIMO channels and elsewhere, and found to be extremely useful. In conducting 
rare-event simulations of systems that involve signal processing operations that 
are mathematically complex, there are two principal issues that contribute to 
simulation time. These have to be dealt with effectively. The first issue concerns 
the rare event itself whose probability is being sought, and this can of course 
be handled by suitable IS techniques. The second concerns the computational 
intensity that accompanies the signal processing. The two are not unrelated. 

In the case of STAP detectors, the chief processing burden is from inversion 
of the sample matrix. It is a daunting task to conduct conventional Monte 
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Figure 1.6: Optimum biasing (scaling) parameters for IS. 



Figure 1.7: Estimated simulation gains. 
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L = 704, N = 352 



Figure 1.8: Simulation gains. 



Figure 1.9: Simulation trials I\ required. 
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Carlo simulations that involve several millions of trials to estimate low FAPs, 
with as many matrix inversions. Assuming that a good IS scheme can reduce 
the number of trials to, say, only a few thousands would still be computationally 
demanding (a case in point being the three thousand 352 x 352 matrices that 
were inverted here). This is the point at which IS departs from conventional 
Monte Carlo in a subtle but important manner. It is almost totally useless to 
run the same random variables through a system in a straight MC simulation. 
With IS however, much can be learnt by repeatedly using the (same) random 
variables. In fact, this is one of the powerful features that adaptive IS (and 
inverse IS) can embed into complex system simulation. 

But how does an IS scheme become effective in the first place? Assume 
that we have a biasing scheme that promises to be effective once the param¬ 
eters of the biasing distributions have been optimized. For large systems (in 
the sense of number of random inputs involved), running truly randomized IS 
algorithms adaptively could become demanding as pointed out above. If sys¬ 
tem performance can be characterized in terms of certain random ‘metrics’ (we 
use the word with a slight abuse of terminology), then these metrics can be 
pre-computed for a given set of input variables, and used repeatedly (which, in 
complex systems such as STAP detectors eases the computational burden) in 
adaptive biasing optimization algorithms. These latter algorithms themselves 
usually require no more than 100 iterations and can be extremely fast. Resulting 
IS simulation gains can be simultaneously estimated and these tell us whether 
we need more or less pre-generated variables to achieve certain accuracies. Ad¬ 
justing this latter number, biasing and system parameter optimization (inverse 
IS) algorithms can be run, once. Thus there is an initial stage of at most a few 
steps during which gains are estimated based on pre-computed metrics and the 
number of these metrics is adjusted. 

All this is not as complicated as it appears. Turning attention to (1.9), 
(1.11), and (1.14), the only two random quantities (or metrics) that are needed 
to estimate the FAP and associated detection threshold are B and D. This is for 
the (/-method. For 2-d biasing the only additional quantity required is the norm 
A of the primary vector, defined just after (1.5) and this adds almost nothing 
to the computation. It turned out that generating K = 3000 random instances 
of B and D was certainly an overkill. If one looks at Figure 1.7, the gain 
provided by IS for estimating a = 10 -6 is about 10 6 . From usual asymptotic 
normality arguments, [5], it follows that about 100 optimally biased trials are 
sufficient to guarantee an absolute estimation accuracy not exceeding 10% with 
95% confidence. For 2-d biasing, the simulation gain will be somewhat lower 
but the essential advantages of the method above remain. That is, handling 
a few hundred inversions (once) is not at all a tall order. This method can 
produce such an avalanche of results that it is tempting to think of it (with a 
slight stretch of imagination) as a ‘turbo-IS’. The above ideas certainly need 
quantification but it is beyond the remit of this short report to delve deeper. 

An interesting observation comes from Figure 1.6, which shows that the 
biasing parameter is very close to unity and has a small spread despite the 
range of FAPs considered. The implication is that the (one-sided) density of 
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the metric D has small variance, presumably owing to the choice of L and 
N. Smaller values of these constants would probably lead to larger spread of 
biasing parameter. In the actual adaptations, a small value of the step-size 
parameter was used to ensure gradual but safe convergence. This explains their 
apparently slow nature as seen in the figures. While configuring results for a 
suite of system parameters, only the first adaptation need be somewhat long; 
subsequent adaptations can be much shorter as they pick up starting or initial 
values from the previous one. 


1.5 Conclusions 

We have made a small inroad into the use of adaptive IS algorithms to char¬ 
acterize a STAP detector. The AMF was used as example and results have 
been very good. The chief reasons for this are that we were able to invoke the 
^-method and inverse IS, find a suitable biasing strategy that could be easily 
optimized adaptively, and find a way around the difficult task of inverting large 
matrices several times (as described above). The hope is that applications to 
other STAP configurations, such as normalized AMF and those that handle non- 
homogenous clutter, will also meet with success. But this remains to be seen as 
we are certainly not in position to predict what subtleties these other detection 
algorithms can throw up. It is clear that IS is still in its infancy, especially 
insofar as use for characterizing modern detection algorithms is concerned. The 
simulation experiments conducted here have raised questions that need to be 
answered subsequently. 


Appendix 

Adaptive algorithms for 2-d biasing The /-function is estimated as 

and its minimization can be carried out using the 2-dimensional adaptive algo¬ 
rithm 

i>m+i = v m - <5J(1.16) 

Here, 5 is a step-size parameter used to control convergence, and m is the index 
of recursion. This is a stochastic Newton recursion. It achieves minimization of 
/ by estimating a solution of 

V/(iz) = (T a I e ) T = 0 

where I a = dl{v)/da and Ig = dl{v)/d0. The estimate of the Jacobian J 
(which is the Hessian matrix of I) is given by 

j _ I ha ho 

\he he 
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where I xy = dl x /dy. It is straightforward to show that the various /-functions 
defined above can be obtained by the notational equations 


4 = E i ,{l{A)WW x } 
Ixx = E+{l(A)WW xx } 
I xy = E ir {l(A)WW xy } 


with various derivatives of the weighting function calculated as 

W a = 

- 


dW 

da 


We = 


dW 

W 


( LN ~j) 


B\W 
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w nn . = 


W e0 = 


W a e = 


d 2 W 
da 2 
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d 2 W 
de 2 

( LN - 2 f)( LN ~ 1 ) + 

d 2 W 

dadO 
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and they can be estimated as in (1.15). The FAP estimator in (1.4) and the 
adaptive biasing algorithm of (1.16) are implemented simultaneously. 
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Chapter 2 


Geometric-Mean and 
Envelope-AMF STAP 
Detectors 


2.1 Introduction 

In the previous chapter we considered the FAP and threshold estimation of AMF 
detector by using IS and inverse IS respectively. Beside FAP, another perfor¬ 
mance measure of a detector is the detection probability, Pd- The higher the 
detection probability for certain SNR, the better the detector. For the AMF de¬ 
tector, the expression for detection probability in homogeneous Gaussian back¬ 
ground is known (see [1]). If the expression for Pd of a detector is unknown 
then it is usually estimated by using standard MC simulation. Analytical ex¬ 
pression for detection probability of AMF in the presence of interfering targets 
is unknown. In this chapter, the performance of AMF for interfering targets 
will be simulated to evaluate its performance. 

As noted in the previous chapter, the AMF detector is similar to the cell¬ 
averaging (CA) detector (see (1.3)). The CA detector take the sums of magni¬ 
tude square of the training data. In the scalar case (only 1 antenna element), 
it has been shown (see Chapter 6 in [5]) that CA detector is better than ge¬ 
ometric mean (GM) detector in homogeneous background, but GM detector 
performs well in the presence of interfering targets. It also has been shown in 
[7] that mean-level CFAR processors including CA preceded by an envelope de¬ 
tector have more robust performance in the presence of Gaussian clutter power 
transitions and interfering targets as compared to the square law detector. 

Following these facts, we propose the GM-STAP detector and envelope-AMF 
detector (E-AMF), array processing versions of GM and envelope CA detectors 
respectively. These STAP detectors should have small loss in homogeneous back¬ 
ground case and will show more robustness compared to AMF in the presence 
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of interfering targets. 

In this chapter, we will give the expressions for GM-STAP and E-AMF de¬ 
tector tests and also prove their CFAR property. In the performance evaluation 
section, FAP and threshold estimation of these detectors by using IS will be 
given, followed by estimation and comparison of the detection probability of the 
proposed detectors with AMF in the homogeneous case and in the presence of 
interfering targets. 


2.2 The GM-STAP detector 

It is well known that the scalar GM detector test can be formulated, in general, 
as 



where U is the output of the test cell (primary data), U {s are the outputs of 
the reference cells (training data), L is the number of reference cells, and Q is 
some multiplier. Following the same formulation, we will introduce GM-STAP 
detector that can be formulated as 



where s^R 1 x and s^R : x(Z) are the outputs of the matched filter with RMB 
beamformer weight R -1 s applied to test cell and reference cells respectively. 
The threshold multiplier for this detector is r] g . The definition of R, s, x, x(l) 
and L are the same as in previous chapter. 

The GM-STAP detector does geometric mean processing of the matched 
filter ouputs of reference cells whereas the AMF detector computes an arithmetic 
average. 

2.2.1 Known covariance 

In the ideal case when the covariance of the data is known, the estim ate R is 
replaced by R in the test of (2.1). Normalizing the test by yields 



where V and {V)} are all independent. This is exactly the same situation as 
for the GM-CFAR detector described in [5]. The FAPs of the GM-STAP and 
GM-CFAR detectors are therefore the same in this case. Furthermore, when the 
covariance is known, the GM-STAP detector has (almost trivially) the CFAR 
property. In any case, such a situation is of no practical interest, for if the 
covariance were known then one would just implement a coherent matched filter 
with fixed threshold. 
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Asymptotic threshold 


It is of interest, however, to determine the threshold of the GM detector as the 
number of training vectors L —> oo. For convenience we consider the square-law 
version of the test in (2.1). Denoting the random part of the (squared) threshold 
by T gm , we have 

r oM w = (n,!, is , ft- 1 ^(oi 2 ) 1/i 

so that 

log t gm ( l ) = \ lo & l st R-~ lx (01 2 

Noting that R —^ R 1 as L —> oo, then s^ R X(Z) stR _ 1 X(Z) in the 
absence of target. Therefore, by the law of large numbers 

log T qm (L) —> Fl{log |s^R _ 1 X(/)| 2 } (in probability) 

i r°° 

= / (log*) ex P -*/ stR " s dx 

s' tt. s j 0 

= - 7 + log(s t R _ 1 s) 

because stR _ 1 X(Z) ~ CA/"i(0, slR _ 1 s) and where 7 is the Euler constant. It 
follows (by Theorem 2.7 of [ 8 ]) that 

T gm (L) —» exp -7 (s^R _ 1 s) (in probability) 

= 0.561459 (s t R~ 1 s) 


- r GM (°o) 

The FAP of the GM-STAP detector in this known covariance matrix case is 
given by 

a GM (L-> 00) = p(|stR- 1 X | 2 >^T GM (oc)|ffo) 

= exp- 0 - 561459 ^ 

Note that the constant 0.561459 is the same asymptotic normalized weight v gm 
obtained for the scalar GM-CFAR detector in Table 6.3 of [5]. Hence, the 
threshold 77 2 required to provide a FAP of 10 -6 is 13.81551/0.561459 = 24.6064, 
in the asymptotic case. 

2.3 The envelope-AMF detector 

The AMF test described by equation (1.3) in the previous chapter involves 
the magnitude square of the matched filter ouputs of primary and training 
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data. We may call this detector as Square Law - AMF (SL-AMF). By a slight 
modification to this detector, taking the magnitude of the filter outputs instead 
of the magnitude squared of the filter outputs, we have the test 

|stR- 1 x||^X)Lls tfi '" 1 x(OI ( 2 - 2 ) 

where the symbols used here have the same meaning as in the previous section 
and Tj e denotes the threshold multiplier for the envelope detector. We call this 
the envelope-law AMF (E-AMF) detector. 


2.4 CFAR property 


In this section an invariance property is established that is of use in constructing 
STAP detection algorithms with FAPs that do not depend on the data covari¬ 
ance R. Although proof of the proposition given here follows the same line of 
argument essentially contained in the exposition of the generalized likelihood 
ratio STAP detector test (Kelly’s GLRT) given in [2], it is outlined here for 
convenience of the reader 1 . Assume, as before, that the primary and training 
data vectors have the same covariance. Consider the variables 

G = s t R- 1 x and G(l) = s t R- 1 x(/) (2.3) 

for l = 1,..., L, that are involved in the AMF detection test of (1.3). Using the 
transformations u = R _1 / 2 s, y = R _1 / 2 x, and y(Z) = R _1 / 2 x(£), leads to 

G = u^R : y and G(l) = u^R _1 y(Z) (2-4) 

where R = R^ 1//2 RR -1 / 2 . The whitened vectors Y and Y (l) are both dis¬ 
tributed CA/"jv( 0,I). It turns out that R has the complex Wishart distribution 2 
CW(L , N; 4-1). Further, a unitary transformation U can be found which rotates 
the new signal vector u into an elementary vector e as 

de = U^u 


such that e = [1,0,... , 0]^ and where d 2 = ||TJdu|| = sfR _1 s. The first column 
of U is the new signal vector u. The remaining columns comprise an orthonor¬ 
mal basis determined, for example, by a Gram-Schmidt procedure. Let z = LPy 


1 II would be helpful for the reader to refer to [2] (see also [1]). We have used several results 
from this now classic paper and have attempted to maintain the same notation. 

2 When X ~ CA'.vfO. R), the (Wishart) matrix W = M(b, XX' = LR has the complex 
Wishart distribution CW(L, N; R) specified by the density 


/w( w ) = i 0 J(R) 


- exp(—tr(R 1 w)),ifRis positive definite 


otherwise 


(2.5) 


where J(R) = n N ( N ~ 1 '>/ 2 fl^Li r U - n + l)(det R)U The covariance estimate R is dis¬ 
tributed as CW(L, N ; -^R). If B is an N x N nonsingular complex matrix, then V = B^WB 
is distributed as CW(L, AT; B^RB). See [4] for more on Wishart distributions. 
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and z(l) = Lby(Z). Applying these to (2.4) yields the variables 


G= 

Jj 


z and G(l) = —e^S 1 z (l) 

1j 


( 2 . 6 ) 


where S = L UtRU. While Z and Z(l) are distributed as CA/"Ar(0,I) and are 
independent, S has the distribution CW(L, N; I). The vectors z and z (l) are 
decomposed as 


and z (l) = 


z aV) 


z b(0. 


where the A components are scalar and B components (N — l)-vector. Corre¬ 
spondingly, S is decomposed as 


5 = V = z(/)z(/) f = 


$AA $ AB 

c c 
^ BA 


(2.7) 


with 


V = S~ 1 = 

The entries of V can be expressed as 


V V 

' aa 1 AB 

V V 

1 BA ' BB 


V 

1 AA 

(^AA ^AB^BB^BA 

V 

1 BA 

= -^BB^BA^AA 

v An 

= v~ x 

' AB 

BA 

V 

' BB 

= +'p- lr p p 

' / AA ; BA ' AB 


as shown in [2] (page 120). Using these definitions and relations in (2.6) gives 

G = — e^5 _1 z 
L 

= £ A A A ~ ^AB^BB Z b') 

d m 

= aaV 


( 2 . 8 ) 


and 


G(l) = — e t 5- 1 z(?) 

Jj 


= T - s AB s-'^ B m) 
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(2.9) 


where 


d 

~L 


VaaV (0 


V = Z A ~ S AB S BB Z B 

y{i)^^ A {i)~s AB s~^ B {i) 


( 2 . 10 ) 


Conditioned on the vectors z B and z B (l), the random variables Y and Y{1 ) in 
(2.10) are (in the absence of target) uncorrelated and Gaussian with zero means 
and variances that can be calculated 3 as 


£ b {|Y| 2 } = l + 

= Z B I 2 ' 11 ) 

using (2.7) in the second step, and 

E B {\Y(l)\ 2 } = l-z B (iyS~ 1 B z B (l), l = l,...,L (2.12) 

with E b denoting conditional expectation. Further, the conditional covariance 
of the variables Y ( l ) is given by 

E B {Y(k)Y(ny} = -z B ( n yS- 1 B z B (k), k^n (2.13) 

Hence the set of conditionally jointly Gaussian zero mean random variables Y 
and {y(Z)| have individual variances and covariances that are functions of the 

random vectors z B and {z B (Z)} . The latter are all jointly independent, each 
being distributed as CAf n-i(0, !)• The probability of any event defined on the 
random variables Y and {Y"(Z)} in (2.10) can thus be determined by perform¬ 
ing an averaging operation over the distributions of z B and {z B (Z)} and this 
probability will be independent of the data covariance R. This statement is also 
true for the random variables G and {G^)]^ in (2.8) and (2.9) with the caveat 
that any constant scaling of these variables should leave the event unchanged. 
The preceding arguments therefore constitute proof of the following 

Proposition 1. Any STAP detection algorithm that uses only the random vari¬ 
ables G and {G^)]^ defined in (2.3) for its description such that the algorithm 
itself is unchanged by arbitrary but equal scaling of all these variables, has a 
FAP which is independent of the target-free data covariance R. 

Proof. As above. □ 

3 See pages 121 and 122 of [2]. 
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The proposition above helps to establish the CFAR property of GM-STAP 
and E-AMF detectors. We can rewrite GM-STAP test in terms of random 
variables G and G(l) by applying (2.3) to (2.1) as 

i g i i % (n, =1 i G (oij (2-i4) 


By applying (2.3) to 
and G(l) 


(2.2), we have an expression for E-AMF in terms of G 



(2.15) 


Scaling both variables G and G(l) with arbitrary but equal scaling factor in 
both above equations will not change the tests. So by using proposition 1, it is 
clear that GM-STAP and E-AMF have the CFAR property. 


2.5 Performance Evaluation 

In this section we present estimates of the false alarm and detection probability 
(Pd) of the proposed detectors for several cases. The FAP of proposed detectors 
cannot be derived analytically so that IS will be used to do the estimation. 
The inverse IS will be executed to estimate the threshold of proposed detector. 
The threshold then will be applied to evaluate the detection probability of the 
detectors. As neither the exact formula nor an approximation for Pd is known, 
it will be estimated using MC simulations. 

2.5.1 False Alarm Probabilities, Thresholds, and Gains 

To estimate the FAP of both proposed detectors, we use ^-method which is 
previously described in chapter 1. The primary data x and training data vectors 
x(Z) are jointly independent and complex Gaussian, sharing a true covariance 
matrix R. Because GM-STAP and E-AMF have CFAR property, we can set 
R = I to simulate the performance of both detectors. In estimating the sample 
covariance matrix for FAP estimation, the training data vectors are assumed 
free of target signal. 

FAP of GM-STAP detector 

To evaluate the FAP of GM-STAP, for the sake of simplicity, we can write the 
test in (2.1) as 4 


BtR-ixl 2 


H i 

H 0 


ILya lx WP 


l/L 


(2.16) 


4 The LHS is, conditioned on the covariance estimate R, an exponential random variable. 
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Using above test, the FAP can be expressed as 


= p |s t R- 1 xr > n, 


2 . 2 
^ <vy" 

9 


IL =1 ^R-'x^i' 


l/L 


= E jP |^|s^R _1 X| 2 > r? g ^ J=1 I^RT'XtOf 
= E{g(X L )} 


H 0 

1 /L 


x L ,// 0 


(2.17) 


By following same derivations as in section 1.3.2, the only changed variable is 
D (see (1.10)). We can redefine D for GM-STAP as 


D g 


(nf =1 istR- i x(oi 2 ) 1/i 

st(R-i)2 s 


(2.18) 


The numerator in the equation above is the right-hand side of GM-STAP test 
replacing the right-hand side of AMF test as in (1.10). Therefore, </(Xl) be¬ 
comes 

S(xi) = e~<° G 

and the FAP expression in term of Dq is 


a = E{e ^ Dg } 


(2.19) 


The weighting function W(X.^; 9), /-function and its derivatives are the 
same as described in the section 1.3.2. The threshold is estimated adaptively 
by using recursive formula described in (1.13). 

Some simulations are run for case L = 128 and N = 64, by applying the 
(/-method. The adaptive IS and inverse IS method are run simultaneously to 
estimate the optimal biasing parameters 9 and thresholds. The scaling is used 
as biasing technique. We simulate FAP from 10 -1 to 10 -7 with logarithmic 
steps of one decade. 

We can see the threshold estimation by using adaptive inverse IS method in 
the Figure 2.1. For each FAP, the threshold reaches its convergence to optimal 
value in 20 iterations. By using these thresholds we get FAP estimates pictured 
in Figure 2.2. The convergence of biasing parameter 9 of each FAP to its optimal 
value can be seen in the Figure 2.3. From this figure we can observe that the 
biasing parameter for each FAP estimate is very close to unity. This means that 
the metric Dq has small variance, probably caused by the choice of L and N. 


Thresholds (GM) 

Figure 2.4 shows the thresholds ry g for corresponding FAP. This can be approx¬ 
imated by a linear interpolation between g g and log 10 a which can be expressed 
mathematically by the relation 

Vg = —l-31og 10 a + 3.2 
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L=128, N=64, Homogeneous Gaussian background 



Figure 2.1: Convergence of rj g . Inverse IS method for GM-STAP detector. 


L=128, N=64, Homogeneous Gaussian background 



Figure 2.2: FAP using IS for GM-STAP detector. 
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L=128, N=64, Homogeneous Gaussian background 



Figure 2.3: Biasing parameter 9 for GM-STAP detector. 


L=128, N=64, Homogeneous Gaussian background 



Figure 2.4: 77 vs. FAP for GM-STAP. 
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The gradient of the linear graph is 1.3 which means if we decrease the FAP by 
geometric factor 10 -1 then the threshold will increase by arithmetic factor 1.3. 
A more accurate relation between r) g and logger can be represented by the cubic 
interpolation 


r]g = -0.011 (log 10 a) 3 - 0.2 (log 10 a) 2 - 2.2 log 10 a + 2 


( 2 . 20 ) 


To be noted is that this interpolation has been obtained for the FAP range 
10" 1 - 10 -7 . 


Simulation gains 

The simulation gains for different FAPs can be seen in Figure 2.5. The obtained 
gains are quite high but lower compared to AMF (L = 128, N = 64) and much 
lower compared to AMF for L = 704 and N = 352. We can expect that for 
higher L and N, the simulation gains of GM-STAP will be higher. However, 
using IS in simulating the FAP of GM-STAP indeed reduces the computational 
load. As an example, for FAP 10~ 6 , the trials needed are only about 2300 while 
MC needs 10 8 . 


L=128, N=64, Homogeneous Gaussian background 



Figure 2.5: Simulation gains T vs FAP for three detectors. 


FAP of E-AMF 

The false alarm probability of E-AMF is also estimated by using IS and g- 
method. To derive the ^-method estimator, for the sake of simplicity, we can 
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express (2.2) as 


IstR^xl 2 I ^(^Etil^'xWl ) 2 

|s t R.- 1 x| 2 f r? e Z 
Ho 

where 

The FAP for E-AMF then can be written as 


a = P^R-'X] 2 > r] 2 e Z\H 0 ) 

= E{P(|s t R~ 1 X| 2 > r) 2 Z\X L , H 0 )} 
= E{g(X L )} 


( 2 . 21 ) 


( 2 . 22 ) 


By carrying out the same steps as in section 1.3.2, we get the expression for 


?(Xl) as 

g{X L )=e~^ DE 

where 

Z 


De = - 


st (R _1 ) 2 s 


The expression of FAP becomes 


a = E{e~^ DE } 

= £ l *{e _T? = De W(X. l - 0)} 


and its estimate is 

Sg= 

The variables needed for adaptive IS like scaling parameter 0,weighting func¬ 
tion W{Xl]9), /-function and its derivatives are the same as described in sec¬ 
tion 1.3.2. To find the optimal biasing parameter 8 we used Newton recursive 
algorithm described in (1.12). The thresholds for each FAP are adaptively esti¬ 
mated by using inverse IS. We may reformulate the adaptive inverse IS algorithm 
described in (1.13) and (1.14) for the E-AMF case as 

Vm-\-l — Vm + 0 V ^ , TYl — 1,2,... (2.23) 

at'iVm) 

where S v is a step-size parameter and the derivative estimator in the denomi¬ 
nator is given by 

a'(Vm) = _^^^e-^W(X L ;0); ~ /* 
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= -2 r/m D e a(rj m ) ( 2 . 24 ) 

with the prime indicating derivative with respect to r/ m . The symbols rj m and 
to correspond to the ?y e and iteration number respectively. 

We test adaptive IS and inverse IS algorithms simultaneously for case L = 
128 and N = 64 to estimate FAPs from 10 -1 to 10”'. The estimated FAPs 
are described in Figure 2.7 by using thresholds found from optimization pro¬ 
cess pictured in Figure 2.6. We can see that the adaptive inverse IS algorithm 
converges in 20 iterations which is quite fast. The biasing parameter (9) opti¬ 
mization process also converges in 20 iterations as pictured in Figure 2.8. As in 
case of GM-STAP detector, the biasing parameters for E-AMF are also close to 
unity and has small spread. 


Thresholds (E-AMF) 

The relation between the thresholds and -logio is considerably more linear in 
this case with gradient 1, as described in Figure 2.9 and expressed as 

Tie = -log 10 a + 2.7 

A more accurate relation is expressed in term of the cubic equation 


r] e = -0.01 (log 10 a) 3 - 0.17 (log 10 a) 2 - 1.9 (log 10 a) + 1.7 


(2.25) 


By using this expression we may have rough estimates of the thresholds for 
lower FAPs. 


Simulation gains 

Figure 2.10 shows the simulation gains by implementing IS simulation with g- 
method. We can observe that simulation gains of E-AMF are slightly higher 
than for GM at high FAPs and almost the same at low FAPs. Compared to 
the AMF, the gains of E-AMF are almost the same at high FAPs but lower at 
low FAPs. As in case of AMF, larger values of L and N may result higher in 
simulation gains. 

2.5.2 P D in homogeneous case 

In this subsection we present the detection probability of GM-STAP and E- 
AMF detectors in homogeneous background and compare them with the AMF. 
We set FAP to 10 -6 , the value of L and N are the same as in FAP simulations. 
All the simulations concerning detection probability are done by using standard 
Monte Carlo (MC) technique, because the probabilities are not very low. The 
number of samples are kept larger than 100/Pd. The detection probability in 
case of non-fluctuating target can be seen in the Figure 2.11. It shows that 
the detection probability of AMF is better than GM-STAP and E-AMF in the 
homogeneous case. For P E equals to 0.5, E-AMF has about 0.18 dB detection 
loss compare to AMF, while GM-STAP suffers 0.3 dB loss. 
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L=128, N=64, Homogeneous Gaussian background 



Figure 2.6: Convergence of rj g . Inverse IS method for E-AMF. 


L=128, N=64, Homogeneous Gaussian background 



Figure 2.7: Convergence of FAP using IS for E-AMF. 
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L=128, N=64, Homogeneous Gaussian background 



recursions 

Figure 2.8: Adaptive IS biasing 6 for E-AMF. 


L=128, N=64, Homogeneous Gaussian background 



Figure 2.9: ?? vs. FAP for E-AMF. 
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L=128, N=64, Homogeneous Gaussian background 



Figure 2.10: Simulation gains T vs FAP for three detectors. 


We also do some simulations to estimate detection probability when the 
Swerling I fluctuating target model is included. To simulate the performance 
because of this fluctuating target model, a combination of ^-method and MC is 
used. This is done to obtain smoother results. The result is presented in Figure 
2.12. We can observe that for Pd equals to 0.5, GM-STAP has loss about 0.3 
dB and E-AMF has loss about 0.15 dB compared to AMF. 

From these simulation results we can conclude that in homogeneous case, 
regardless the fluctuation of the target, both proposed detectors have small 
detection loss compare to AMF. 

2.5.3 Pd in the presence of nonhomogenities 

In this subsection we will show the simulation results that compare the perfor¬ 
mances of GM-STAP and E-AMF with AMF in the presence nonhomogenities. 
The nonhomegenity we consider here is the presence of interfering targets. The 
power of the interfering targets are taken to be the same as the power of the ac¬ 
tual target 5 . Figure 2.13 shows the comparison of the performance of GM-STAP 
and AMF in the presence of 2 interferers and the target signal is constant,i.e. 
non-fluctuating. At Pd = 0.5, AMF has detection loss about 0.98 dB w.r.t 

5 It is assumed that interferers have the same angle-Doppler properties as those of the 
primary target. 
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L=128, N=64, Swerling 0 target, 2 interferers 



Figure 2.13: P D of AMF, GM-STAP, and E-AMF. 


GM-STAP and about 0.73 dB compared to E-AMF. This detection loss will 
increase for higher Pd- For example at Pd = 0.9, AMF suffers loss about 1.9 
dB w.r.t GM and 1.6 dB w.r.t E-AMF. From the graph, we can also see that the 
detection loss of E-AMF compared to GM is quite constant i.e about 0.25-0.3 
dB. 

For case 3 interferers, the losses of AMF are even higher. At Pd = 0.5, AMF 
has loss about 2.5 dB w.r.t GM and 1.9 dB w.r.t to E-AMF. For an extreme 
case, i.e. at Pd = 0.9, AMF has loss about 9 dB and 8 dB compared to GM 
and E-AMF respectively. The E-AMF itself has loss 0.58 dB (at Pd = 0.5) and 
0.85 (at Pd = 0.9) dB compared to GM. 

These facts show the robustness of GM and E-AMF detectors compared to 
AMF. In the case of Swerling I fluctuating target model, both detectors still 
show their robustness compared to AMF as can be seen in the Figure 2.15. We 
can conclude that in the presence of interfering targets, GM-STAP performs the 
best while AMF is the worst. The E-AMF detector also shows its robutness in 
these simulations and has small detection loss compared to GM detector. 


2.6 Conclusion 

In this chapter we introduced the GM-STAP and E-AMF detectors. We show 
that IS methods work well in estimating the FAP and determining the thresholds 
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Table 2.1: Comparison of three STAP detectors. 


Q 

AMF-S 

AMF-E 

GM-STAP 

Is it CFAR? 

Yes 

Yes 

Yes 

homogeneous Pd 

- 

very small loss 

very small loss 

Robust to interferers? 

No 

somewhat 

most 


of GM-STAP and E-AMF detectors. However, IS simulation gains for GM- 
STAP and E-AMF are not very high, and they can be increased by developing 
better biasing methods which is subject for future research. 

Simulation results show that the GM-STAP and E-AMF detectors have quite 
small loss in detection probability compared to AMF, but they have very good 
performance and show robustness in the presence of interfering targets. The 
AMF fails to maintain its robustness in this situation. These findings are sum¬ 
marized in Table 2.1. 
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Chapter 3 

Further results on IS for 
STAP detection 


This third report (or chapter) begins with some additional results on the E-AMF 
and GM-STAP detectors for the case of L = 704 and N = 352. These complete 
our investigations into the two detectors and complement parallel results that 
were obtained previously for the square-law AMF detector with the above L 
and N values. Comparisons with the latter detector have also been made. 

The results are contained in Figures 3.1 - 3.5 for the E-AMF detector. Inverse 
IS results are in Fig. 3.1 and 3.2. Estimated optimum biasing parameters and 
the /-functions are in Fig. 3.3 and 3.4 respectively. IS gains estimated during 
simulation are shown in Fig. 3.5. Corresponding results for the GM-STAP 
detector are in Figures 3.6 - 3.8. Thresholds for the 2 detectors (including 
AMF) as functions of FAPs are in Fig 3.9. The cubic interpolations 

g E = 0.01a; 3 - 0.18647a; 2 + 1.8822a: + 1.725467 

and 

rj G = 0.011813a; 3 - 0.21985a; 2 + 2.22917a; + 2.040834 

where x = — log 10 a can be used to determine thresholds for the E-AMF and 
GM detectors in the FAP range 10 _1 - 10~' respectively. Estimated IS gains as 
functions of FAPs are shown in Fig 3.10. It is noted that simulation gains for 
the 3 detectors are of the same order. 


3.1 General remarks on IS 

Point of application of biasing 

Some comments of a general nature regarding application of IS to signal process¬ 
ing algorithms are made here. They give some insight into the thinking behind 
our efforts to simulate the above mentioned algorithms. The simulation proce¬ 
dures described thus far in the previous chapters are applicable to any detector. 
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Envelope Detector with N = 352, L = 704, K=5000 



Figure 3.1: Inverse IS thresholds. 
Envelope Detector with N = 352, L=704, K=5000 



Figure 3.2: FAP estimates. 
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Figure 3.5: IS gain estimates. 

Geometric Mean detector, N = 352, L = 704, K = 5000 
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Figure 3.6: Inverse IS thresholds. 





















































Geometric Mean detector, N = 352, L = 704, K = 5000 



Figure 3.7: FAP estimates. 

Geometric Mean detector, N = 352, L = 704, K = 5000 



Figure 3.8: IS gain estimates. 








































N = 352, L = 704, K = 5000 
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Figure 3.9: Thresholds as functions of FAP. 
N = 352, L = 704, K = 5000 



Figure 3.10: Gain comparisons as functions of FAP. 































The IS biasing is performed on the input random variables and we refer to this 
as input biasing. Accuracy of the IS estimates and resulting simulation gains 
will of course depend on the particular detection algorithm under study. It is 
generally true (and intuitive) that better estimator performance can be obtained 
if IS biasing can be carried out closer to the point in the processing chain of the 
detector where the actual (rare event) decisionmaking is done. This of course 
necessitates knowledge of density functions of the processes at the point where 
biasing is to be implemented. Often, input stochastic variables may have un¬ 
dergone transformations whose results are difficult to characterize statistically 
in analytical form, and we have to rely on the general method above. However, 
when such transformations can be characterized, then IS should be carried out 
using the modified processes. Therefore, it is desirable to perform biasing and 
IS as close to the final decisionmaking point as permitted by availability of 
knowledge of probability density functions. 

Another approach to biasing that is sometimes possible is to perform a series 
of (linear and/or nonlinear) transformations of the input processes as if one were 
carrying out a mathematical analysis of the algorithm. The transformations are 
carried out until the point beyond which it may not be possible to determine 
the density functions of the transformed processes without considerable mathe¬ 
matical effort. Biasing is then performed at this stage, the hope being that the 
g-method becomes applicable. This procedure may produce higher simulation 
gains than simple input biasing which, of course, is the easiest to implement. It 
is illustrated in some of the following examples. 


3.2 The NMF STAP detector 


The use of IS to characterize the normalized matched filter (NMF) detector and 
some of its variants is discussed in the rest of this chapter. Expressions for the 
FAP and detection probability of the NMF detector operating in homogeneous 
Gaussian clutter are known in closed form. However, much can be learnt about 
biasing for this class of detectors from the exercise of applying IS to estimate its 
FAP. The variants are the normalized adaptive matched filter (NAMF) and the 
low-rank NAMF detectors. All these detection algorithms have been treated in 
detail in [3], [15], and other papers. 

The NMF detection test is given by 


A 


NMF = 


|stR 1 x | 2 
(sfR _ 1 s)(xtR _ 1 x) 


H 1 

5 : 

H 0 


(3.1) 


following the usual notations. Its FAP in Gaussian interference is known, [15], 
and particularly easy to derive. It differs from the matched filter (for known 
interference covariance matrix R) in the normalization term which is the second 
one in the denominator. The FAP of the detector is given by 


“nmf 
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and it has the CFAR property of being invariant to the interference covariance 
matrix R. Furthermore, as is evident from (3.1), the FAP is also invariant to 
any scaling of the primary and secondary data. 

3.2.1 The FAP of the NMF detector (using IS) 

We begin with a simple (re-) derivation of the FAP of the detector. It helps to 
illustrate an elegant aspect of IS, by which it is sometimes possible to derive a 
perfect estimate of a rare-event probability. Furthermore, although this problem 
can be analyzed completely, it is an application of the procedure described in 
the last paragraph of the introduction. 

It is assumed that x ~ CMn{ 0, R). With a whitened data vector defined as 
xi = R^ x / 2 x and a transformed steering vector Si = R _ 1 / 2 s, the test statistic 
of (3.1) takes the form 


Anmf 



x|s! s|x! 1 
ll S l|| ll S l|| xjxi 


xjuiujxi 

xjx! 

y*yi 

xjxi 

N 2 

llxill 2 


(3.2) 


where Ux = si/||si || is an IV-dimensional unit vector and yi = ujx! a scalar ran¬ 
dom variable. Using for example a Gram-Schmidt procedure, an IV-dimensional 
basis can be formed by determining N — 1 other unit vectors in the orthogo¬ 
nal subspace of ui. Denoting the former by U;, i = 2we define the 
corresponding random variables y t = ujxi. Then 


N N 

\yi\ 2 + J2 i= Jy^ 2 = 


i=l 


N 


= X 


55 u * u ! x i 


i =1 

= xjlxi 
2 


= xi 
= l|xi 


(3.3) 
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with the {y;}^ being i.i.cl. and distributed as CA/"i(0,1). 
becomes 


AnMF 


M 2 

ll/ll 2 + E<^2 \Vi\ 2 


H 1 
Ho 


V 


The test therefore 
(3.4) 


which can be put in the form 




H 1 

- ^ 
Vi Ho 


Vo = 


l-rj 


(3.5) 


where 

u=\yi\ 2 and = \yi\ 2 , i = 2,...,N 
The FAP of the NMF detector is then 


ry 

U NMF 


P(u>io'Z i=2 V i ) 

e{p(U > r,oJ2* =2 V * |v 2 ,...,Viv)} 


(3.6) 


The first line above has exactly the same form as the FAP of a CA-CFAR 
detector as U and each V) are i.i.cl. (unit) exponential random variables; the 
formula for this probability is well known. To re-derive the latter, a g-method 
estimator can be written as 


K 


NMF 


= -W(V 2 ,...,V N ) 


Scaling each V) with a and using the resulting weighting function 


W(v 2 ,.. .,v N ) = a 




yields the estimator 


ry 

U NMF 


1 

K 


K 

Y a N ~ X e - (,? ° + 1 - 1 /a) S~ Vi 

1 


(3.7) 


(3.8) 


(3.9) 


If the scaling factor is chosen as a = 1/(1 + rj 0 ) in this estimator, then 


ry 

U NMF 



(1 — rj) N a constant 

ry 

U NMF 


(3.10) 


the last step following from the fact that the variance of this unbiased estimator 
is zero. 


44 



3.2.2 FAP estimation by IS: rotation of primary vector 

If we want to estimate « NMF using input biasing, then it is clear from the 
test statistic in (3.1) that a simple scaling of the elements of the primary data 
vector x will be useless. A form of biasing can however be developed using 
the (well known) fact that the test statistic is actually a cosine-squared one, 
the concerned angle being that between the transformed steering (or unit) and 
whitened primary data vectors Si (or Uj) and xj respectively, as evident from 
(3.2). The frequency of false alarm events in a simulation can therefore be 
increased by biasing xi so as to decrease the angle toward zero or increase it 
toward n. This can be accomplished by a rotation of Xi which, since it is being 
assumed for this detector that the actual data covariance matrix R is completely 
known, is equivalent to input biasing (of the primary data vector x). 

Consider biasing the whitened primary data vector xi by rotation with 
an N x N matrix A. The biased vector is Axi and will be distributed as 
CA/"at(0,R*) with covariance matrix R* given by R* = AA^, since xi ~ 
CA/"/v(0, 1 ). The weighting function is then given by 


W(xi;A) 


/(* i) 

/*( x i) 

|R*I exp ( - x{(I - R* 1 )x i ) 


(3.11) 


where | • | denotes matrix determinant. The FAP estimator for the test in (3.2) 
then takes the form 


K 


a 


NMF 


= ^T! 1 (l u i Xi l 2 > »? ll Xi ll 2 ) W(*i; A); ~/* (3.12) 


The problem now centers around determining an effective rotation or biasing 
matrix A. There is clearly an unbounded number of rotation matrices to choose 
from. To narrow the choice, it is most convenient to make the matrix dependent 
on a single (biasing) parameter, say a, and search for the optimum value of the 
latter by minimizing the associated /-function which is given by 

1(a) = £*{l(A) 1T 2 (x i; A)} 

= £{l(A)ir(x i; A)} (3.13) 

where A = l(|u{xi| 2 > 77 11 x 1 11 2 ) is the false alarm event. In actual simulation 
we will determine a minimizer for an estimate of 1(a), that is 


1 

a aP t = argnun (f(a) = — ^ l(A) W 2 (x.\\ A); ~ /*) (3.14) 

1 

The form or structure of A has to be now decided. We assume that A = I + aT 
for 0 < a < 1 where T is some N x N matrix. The biasing matrix A is to be 
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constructed so as to introduce correlation amongst the components of Xi in a 
controlled manner. To understand how this can be achieved, consider first the 
simple case where the unit vector U! is given by 


Ul 


1 

Vn 


I 1 \ 
V 1 


This can happen if R = I and the original steering vector s has all equal 
elements. Suppose (based on heuristic reasons) the matrix T is set to 


Then 


(° 1 

1 0 


T = 


n 


v i 


o/ 


A = 


l 1 


a 
a 1 


a \ 

a 



(3.15) 


It is noted that the above matrix does not provide a norm-preserving rotation. 
We avoid simulating with a = 1 since this obviously leads to a singular matrix. 
When a = 1, the test statistic becomes unity and all the elements of the biased 
vector Ax! are equal to JJiLi x u where aq.j’s denote the elements of the unbiased 
x-| . That is, the fully biased vector would be aligned in the direction of Ui or 
opposite to it, in C N . When a = 0, no biasing or rotation takes place. Estimating 
an optimum value of a that maximizes the simulation gain for given threshold 
77 is an easy problem. 

For the case R/Iwe consider a unit vector U! having elements denoted by 
ti, i = 1,..., N. By analogy we set 


T = 


/ *i — l 
t2 


1 1 

* 2-1 




\ *iV • *AT 1 / 


(3.16) 


Again, when a = 1 the biased vector is just ( ’Yhi —1 ^li ) 111 anc ^ collinear with 

Ui. 

An interesting phenomenon takes place with this biasing scheme. The weight¬ 
ing function in (3.11) depends on the matrix A which in turn depends (apart 
from on the biasing parameter a) on the elements of the unit vector through 
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NMF detector, N= 128, K = 5000 



Figure 3.11: Optimum biasing parameter for the T matrix in rotation biasing, 

R = I. 


the T matrix. The IS estimator in (3.12) estimates a NMF = i?{l(A)} which, 
from (3.10), is known to be independent of the particular unit vector Ui be¬ 
ing used. Hence, if the estimator S NMF is a reasonably good one, then it will 
be largely unaffected by the choice of the biasing matrix A (and hence of T). 
However, the gain of the proposed IS scheme depends on the /-function given in 
(3.13) and this depends on the weighting function and therefore on the matrix 
T. Thus the maximum gain for the optimized IS scheme will be affected by 
the choice of unit vector u x . Since Ui = R _ 1 /, 2 s/||R _ 1 / 2 s|| it follows that the 
IS performance of this rotation biasing scheme depends on the data covariance 
matrix R in force; this is despite the fact that the detector FAP is independent 
of R. Decoupling the biasing matrix A from Ui is of course not possible. Thus 
there may exist some R which gives a best IS simulation gain for given detector 
constants rj and N. Furthermore, it may be possible to achieve some invariance 
of IS gain to covariance matrix R by using a matrix normalization for the T 
or A matrices. These issues are not investigated here and in the interests of 
expediency we present results obtained by simulating the R = I case. We will 
visit rotation biasing later in this investigation. 

The simulation results are in Figures 3.11 and 3.12 which show the estimated 
optimum biasing parameter and IS gains as functions of FAP respectively. As is 
evident from Fig 3.12 the simulation gain is not very high, being approximately 
4 x 10 4 at a FAP of 10" 6 . 
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NMF detector, N = 128, K = 5000 



Figure 3.12: Simulation gain for rotation biasing, R = I. 


3.3 The NAMF STAP detector 


This section deals with the fast estimation of FAP for the NAMF detector. 
Developing an IS procedure for this is made easier by first setting up an alternate 
IS mechanization for FAP estimation of the AMF detector studied in Chapter 1. 


3.3.1 FAP estimation for the AMF revisited 


The alternate procedure for the AMF detector uses the g-method with two- 
dimensional biasing. 

Reproducing (2.8) and (2.9) on page 20 we have 


and 


where 


G(l)=^-P AA y(l) 

y = Z A- S AB S BB Z B 
y(l)=z A (l)-S AB S-' B z B (l) 


(3.17) 


(3.18) 


(3.19) 
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The AMF test of (1.3) on page 3 can therefore be written as 


\y \ 2 


Hi 

H 0 


?Ei»(oi s 


i=i 


(3.20) 


We know from Section 2.4 that, conditioned on z B and {z B (Z)} (or, for short, 

the B-vectors), Y and {Y”(Z)} are zero mean uncorrelated Gaussian random 
variables with variances and covariance as in (2.11), (2.12), and (2.13). From 
[2] it is known that 

L L-N +1 

EmoM e moi 2 (3-21) 

i=i i=i 

where the w(l) are i.i.d. each with distribution CA/"i(0,1). Moreover, Y/M^ 2 
is, conditioned on the B-vectors, also distributed as CA/”i(0,1) where 


m b = i + v b 

and 

E b =z s(Eti Z H^ z B( / ) t ) Z B (3-22) 

The test in (3.20) then takes the form 

H L-N +1 

u ^ v' E (3.23) 

H 0 r-f 


where U = \Y\ 2 /M B and {C/(Z) = |u;(/)| 2 }^ JV+ 1 are all unit exponential and 
i.i.d., and rj = rj/(LM B ). Once again, this is in the form of the usual CA- 
CFAR test when it is conditioned on the B-vectors. Hence the FAP of the AMF 
detector can be written as 


where 


'AMF = ■ P ( C, 2’>'EJ 1 f'fO) 

—»L—iV+1 

U >ri'Y , l=1 

= E{g( S B )} 

s(s B )= 1 


B-vectors 


)} 


(3.24) 


(3.25) 


[1 + v /(lm b )] l - n + 1 

We can therefore estimate the FAP using the ^-method with an IS simulation 
that biases the B-vectors. This estimator is 


q amf 


1 

K 


E<7(£ s )IP( z s , z b J; 


1 


~/* 


(3.26) 
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AMF detector, N = 64, L = 128, K = 10000 



Recursion Index n 


Figure 3.13: Optimum scaling for primary B-vectors. 

where z bl = 

Biasing of the B-vectors contained in S B must produce an increase in the 
value of the (/-function in (3.25). This means that S B must be made to increase, 
which can easily be accomplished by scaling up the primary B-vector z B and 
scaling down the secondary B-vectors z fi ( l ). A 2-d biasing scheme results, which 
needs to be optimized adaptively. As in Section 1.3.1 on page 3, the primary 
and secondary scaling parameters are chosen as a 1 / 2 and 0 1 / 2 respectively, where 
a > 1 and 0 < 9 < 1. The weighting function is easily shown to be 

W ( z b’ z bl) = a N ~ 1 e L(N ~ l) exp ( - z^z B (l - 1/a)) 

• exp ( - (1 - 1/9) z b (0 tz s ( Z )) ( 3 - 27 ) 

The rest of the optimization procedure is as discussed in Section 1.3. 

The results obtained from this 2-d IS simulation have been compared with 
the gf-method used with input biasing, the technique which was developed in 
Chapter 1. Figures 3.13 - 3.16 contain these and are self-explanatory. From 
Fig 3.16 it is clear that applying the (/-method with biasing of B-vectors is a 
more powerful IS scheme than the corresponding method used with input biasing 
of all secondary vectors, which was the subject of Chapter 1. This illustrates 
the points made in Section 3.1. 
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AMF detector, N = 64, L = 128, K = 10000 



Figure 3.14: Optimum scaling for secondary B-vectors. 
AMF Detector, L = 128, N = 64, a g = 10~ 6 



Figure 3.15: IS gain surface for 2-d biasing of B-vectors. 
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AMF detector, N = 64, L = 128, K = 10000 



Figure 3.16: IS gain comparison for both g- methods. 


3.3.2 IS for the NAMF detector 

The NAMF detection statistic is given by 

Is^R^xl 2 1 


(stR - 1 s) (xtR _ 1 x) h, 
Using (3.23), this can be written as 

u tb rj 


H i 


L-N+l < XjM 

E «(/) Hp LM 

1=1 


(xtR-^) 


(3.28) 


(3.29) 


B 


With the transformations in Section 2.4 on page 19, the normalization term in 
the RHS of the above equation becomes 


x^R 1 x = Lz^S 


(3.30) 


with S defined in (2.7) on page 20. From the discussion on pages 120 and 
121 of [2], it turns out that the quantity z Jf S~ 1 z, denoted therein as S, can be 
expressed as 




-l 


\y\ 2 


E W )\ 2 

1=1 
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\y \ 2 


L-N+l B 

E KOI 2 

1=1 


= M 


u 


B L-N+l B 

E «(0 

Z=1 


(3.31) 


using the definitions of Section 3.3.1. Combining (3.29), (3.30), and (3.31) yields 
the test in the form 


H i 

u ^ ?? 0 
H 0 


B 

l + s 


B 


L-N+l 

u (o 


i-i 


(3.32) 


where r] 0 = rj/(l — ij). Conditioned on the B-vectors this is again a CA-CFAR 
test resulting in the FAP 


“namf = e {9 n (Z b )} ( 3 - 33 ) 

where 

3 n ^ b )= [1 + 7/ 0 E b /(1 + S B )] i_iV +i (3 - 34) 

The FAP estimator is then 

1 K 

“namf = ^( z b> z bz,) ; ~ f* (3.35) 

i 

with the weighting function W being the same as in (3.27) excepting that the 
roles of the two biasing parameters are exchanged, that is, 0 < a < 1 and 9 > 1. 
Simulation results are shown in Figures 3.17 - 3.22. IS results for threshold 
estimation have been compared with the analytical expression available in [15]. 
A quick way of determining the threshold for the NAMF detector is using the 
cubic interpolation polynomial : 

%amf = 7 - 96 x iO -5 ^ 3 - 0.003116 a; 2 + 0.0723 x + 7.77 x 10“ 5 (3.36) 

where x = — log 10 a. 


3.4 Conclusion 

In this chapter we have introduced and developed IS techniques suitable for 
analysis of the NMF class of STAP detectors. In particular, the method based 
on rotation of the primary data vector could be an effective solution for the 
difficult problem of characterizing the low-rank version of NAMF detectors. 
This is a matter for further investigation. Two dimensional biasing methods 
have been used to produce new IS results for the AMF and NAMF algorithms. 
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Estimated a Threshold estimates t| 


NAMF detector, N = 64, L = 128, K = 10000 



0 50 100 150 

Recursion Index n 


Figure 3.17: Thresholds through inverse IS. 


NAMF detector, N = 64, L = 128, K= 10000 



Figure 3.18: Optimum scaling for primary B-vectors 





















NAMF detector, N = 64, L =128, K = 10000 



Figure 3.19: Optimum scaling for secondary B-vectors. 
NAMF detector, L = 128, N = 64, a = 1CT 6 

x 10 5 



Figure 3.20: IS gain surface for 2-d biasing of B-vectors. 
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NAMF detector, N = 64, L = 128, K= 10000 



Figure 3.21: Comparison of simulation and numerical integration for thresholds. 
NAMF detector, N = 64, L = 128, K= 10000 



Figure 3.22: IS gain as function of FAP. 
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Chapter 4 


Geometric mean and 
envelope-law NAMF STAP 
detectors 


4.1 Introduction 

Following the developments of Chapter 2, we propose two NAMF detector vari¬ 
ants. These are the envelope and geometric-mean based versions. They are 
referred to hereafter as the E-NAMF detector and the GM-NAMF detector. 
In this chapter, their CFAR property is established and threshold settings for 
the detectors for specified false alarm probability (FAP) is accomplished using 
fast simulation based on importance sampling. The two variants show robust¬ 
ness in detection probability in the presence of interfering targets contaminating 
the training data while having zero loss in homogenous Gaussian interference 
compared to the usual square-law NAMF detector. This is shown for both fluc¬ 
tuating and non-fluctuating target models. Performance comparisons with the 
AMF variants are provided. 


4.2 E-NAMF detector 


This detector is constructed in the same way by which we defined the E-AMF 
detector. Consider the test 


| S tR -1 x| 

\J x’R- 1 x H o 


|Ei st * lx (oi 

^ i=i 


(4.1) 


which we refer to as the E-NAMF detector. As before, we expect this detector 
form to display some robustness against interfering targets in the secondary 
data compared to the usual NAMF detector. 
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4.2.1 Asymptotics 


The asymptotic FAP behavior of this detector as L —> oo (or known covariance 
matrix R) is obtained from 


a. 




= P 


|sfR~ 1 X | 2 
XtR IX 



- (M*) 


N-l 


(4.2) 


the first line above following from convergence arguments and the law of large 
numbers applied to the independent and identically distributed sequence {s^R _1 X(Z)}, 
the second by noting that sTR _1 X(Z) is distributed as CM i(0, s^R^s), and the 
third from the fact that the second line represents the FAP (1 — v) N ~ l of a nor¬ 
malized matched filter (NMF) STAP detector (the statistic of which is shown 
in (3.1)) with threshold v ([15]) where v = 77 ^ 7 r/4. The asymptotic threshold 
can be calculated from the above expression. 


4.2.2 Alternate form and CFAR property 


Using (3.18) and (3.19) of section 3.3.1, the E-NAMF test of (4.1) can be rewrit¬ 
ten as 


\y\ 

V x'R~ ! x H o 


iEw‘)i 


1=1 


(4.3) 


Then, combining (3.30) and the first line of (3.31) and substituting for x^R : x 
in the above yields the test 


I l!>l I |!>(i)l (4.4) 

\J le b + L l»l 2 /Ef.i l#(()l 2 l ~' 

So if we replace y and y(l) by G and G(l) as defined in [10] respectively in the 
above, then the test remains unchanged. This implies, from the proposition of 
[10] and the arguments leading up to it, that the test is CFAR. 


4.3 GM-NAMF detector 


The geometric-mean version, referred to as the GM-NAMF detector, is defined 
as 


|stR- 


.r 


Hi 

^ V g 


IT, , l st R- lx (0l 


1/L 


(4.5) 


where is the threshold. The square-law version is identical, being just a 
square of the above test expression. 
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4.3.1 Asymptotics 

The asymptotic FAP for this detector is given by 


= P( log 


IstR^XI 2 

XtRA 


> 2l0g??g 


^{loglstR-ixa)! 2 } 


, IstR^XI 2 

= P [ 1<>8 X1R--X a21 ° S ”» 


- 7 + logfs^R ^t) 


IstR-iXI 2 o ^ , , 

= PI XtR-.X >V 7s 'lVs 

N—l 


= (l-V 2 9 e~-) 
where 7 is the Euler-Mascheroni constant. 


(4.6) 


4.3.2 Alternate form and CFAR property 

Just as for the E-NAMF detector, the GM-NAMF test of (4.5) can be rewritten 


as 


and then as 


\y\ 


V xtR 


\y\ 


H 1 

^ Vg 


! x H 0 


nii/(oi 


1=1 


l/L 


ylLZ B +L\y\*/'£tL 1 \y(l)\ 


H 

2 H ° 


^ Vg 


ni»(oi 


1=1 


l/L 


(4.7) 


(4.8) 


For the same reasons as for the E-NAMF detector, the GM-NAMF test is also 
CFAR. 


4.4 FAP estimation 

With some straightforward manipulations, the E-NAMF test in (4.4) can be 
rewritten as 


\y\ 


where 


,2 5 

1 Ho LC e 

(It, MO O' 

(4.9) 

n 2 

C e = 1 -A 

(£f=*|y(OI) 2 

(4.10) 

L 

Ehmi 2 
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and the GM-NAMF in (4.5) as 


\y\ 


7 2 

S lg 

Ho 


2 ^ 


Co 




l/L 


where 


2 v n*=i \y ( l )\' 2 


Cg = 1 - L77 , 

5 Ef=i |y(OI 


(4.11) 


(4.12) 


In writing (4.9) and (4.11) we have assumed that C e > 0 and C g > 0. This 
will be guaranteed if r] e < 1 and y g < 1, respectively, and can be seen by 
an application of Holder’s inequality to the sums in (4.10) and using the fact 
that geometric means are smaller than arithmetic means in (4.12). As shown by 
simulation results, a wide range of FAPs is achieved for values of both thresholds 
less than unity and therefore we assume this restriction to hold. For brevity we 
denote the RHSs of both (4.9) and (4.11) by S^T, where D is a function of 
( y(l) • • • y(L) ) = y|) and corresponding constants that depend on the 
detector in question. With this (generic) notation, the FAP can be written as 

a = P(\Y\ 2 >Z b D) 

= E B {P(\Y\ 2 >E B D\B-vec)} 

= E B {E{P(\Y\ 2 >Y B D\B-vec,Y L )}} 

= E B {E{e~^B D ^ 1+ ^B^ |B-vec)}} 


where 


= E B {E{g(Y B1 Y L ) B-vec}} 


fl (E R ,Y i ) = e- E B D/( 1 +E B) 


(4.13) 


(4.14) 


and E b denotes expectation over the distribution of the B vectors and E the 
expectation over the conditional distribution of Y^. The fourth line above 
follows because conditioned on the B vectors, Y and Yl are independent and 
|Y| 2 /(1 + H s ) is unit exponential. 

From (2.7) on page 20, y and y(l ) are 


y= Z A-Y, l= 1 Z A( l ) Z B ( l ^ S B 1 B Z 

y(l) = z A (l) - V L =1 z A {i)z B {rf S bb „ b 


B 

1 


(4.15) 


MO 


Then the FAP in (4.13) can be further written as 
a = E B {E{g(Y B , Y L )| B-vec}} 
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g(S B , Y l ) f(Y L B-vec) /(B-vec) dY L dB-vec 


= j j g(E B ,Y L )f(Y L , B-vec) dY L dB-vec 

= E{g(E B ,Y L )} 

= E{g(X B ,Y(l),...,Y(L))} 

= E{gp B ,Z A (l),...,Z A (L))} 

= E{g(Z B ,Z AL )} 


I / 9(s - 


z al) /(zal, B-vec) dz AL dB-vec 


J J g(E b , z AL ) /(zal) /(B-vec) dz A L dB-vec 
j J g{ E b ,zal) f*(z AL ) f (B-vec) dz AL dB-vec 

J J g{ S B , zal) W(zal) /*(zal) /(B-vec) dz A ^ dB-vec 


= E*{g(Y B1 Z AL )W(Z AL )} (4.16) 

where g{E B > Zal) = g(Y B ,Y] j ) and is given by (4.14). We perform IS only on 
the variables in Zal = (-Za(I) • • • xTa(T)) t , for simplicity. Noting that Zal ~ 
CAf l( 0,1), scaling down Zal with parameter 9 leads to the weighting function 

W{z AL ) = 0 L e~Ei z f A {i)z A {i)(i-i/8) (4.1 7 ) 


The IS estimator for FAP is therefore 
1 K 

a = — Zal) W(zal)] W ; Zal ~ /*(z A l), B-vec ~ /(B-vec) 

(4.18) 


4.5 Simulation results 

The adaptive IS implementation is as usual. Simulation results for the FAP 
range 10“ 1 — 10 -7 are shown in Figures 4.1 - 4.14, for both ENAMF and 
GNAMF detectors. The parameters used in the simulations are L = 128, N = 
64 and the number of IS trials is K = 50000. Moreover, for L = 128, N = 64, 
thresholds for these detectors can be determined using the cubic interpolation 

rj e = 9.704 x 10 -4 # 3 - 1.85 x 10 " V + 0.1595a; + 0.15613 (4.19) 


61 



(4.20) 


rig = 1.1439 x lirV - 2.148 x 10 " V + 0.188a; + 0.184 
where x = — log 10 a. 


E-NAM F Detector, N = 64, L = 128, K = 50000 



Figure 4.1: Inverse g-method thresholds. 


4.5.1 Pd in homogeneous case 

In this subsection we present the detection probability of the E-NAMF and GM- 
NAMF detectors in homogeneous background and compare it with that of the 
NAMF detector. We set FAP to 10 -6 , the value of L and N are the same as in 
FAP simulations. All the simulations concerning detection probability are done 
using Monte Carlo (MC) technique and the number of trials used is 100000. The 
detection probability in case of non-fluctuating target can be seen in Figure 4.15. 
In this figure is also plotted the detection probability for the NAMF detector 
obtained using the analytical formula of [15]. The figure shows that the three 
detectors perform almost the same in the presence of homogeneous clutter. In 
fact, the E-NAMF and NAMF perform exactly the same while the GM-NAMF 
has a very small loss compared to the other two detectors. Simulations are also 
done to estimate detection probability when the Swerling I fluctuating target 
model is used. The result is presented in Figure 4.16. Similarly, in this case 
the detectors have almost the same performance. From these simulation results 
we can conclude that, in homogeneous case, regardless of the target fluctuation, 
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Figure 4.2: FAP estimates. 



0 50 100 150 

Recursion Index n 


Figure 4.3: Gain estimates. 
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Figure 4.4: I g function estimates. 

E-NAM F Detector, N = 64, L = 128, K = 50000 



Figure 4.5: Optimum biasing parameter estimates. 
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E-NAM F Detector, N = 64, L = 128, K = 50000 



Figure 4.6: Threshold as function of FAP. 


E-NAM F Detector, N = 64, L = 128, K = 50000 



Figure 4.7: Gain as function of FAP. 
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GM-NAMF Detector, N = 64, L = 128, K = 50000 



Recursion Index n 

Figure 4.8: Inverse g-method thresholds. 


both the E-NAMF detector and the GM-NAMF detector perform as well as the 
NAMF detector. 

4.5.2 Pd in the presence of nonhomogeneities 

In this subsection we will show the simulation results that compare the perfor¬ 
mance of the envelope, geometric mean and square law NAMF detectors in the 
presence of nonhomogeneities. The nonhomogeneities considered here are two 
interfering targets in the secondary data. In the simulations we assumed the 
interfering targets to have the same power and the same steering vector as the 
actual target. Figure 4.17 shows the comparison of the performance of E-NAMF, 
GM-NAMF and NAMF in the presence of 2 interferers when the target signal 
is constant, i.e. non-fluctuating. In Figure 4.18 is plotted the detection proba¬ 
bility when both the target and the interferers are assumed to be fluctuating. 
For both target models, it can be seen that the NAMF detector, for a fix Pd, 
has a significant detection loss compared to the envelope and geometric mean 
variants and that the GM-NAMF performs slightly better than the E-NAMF 
detector. For example, at Pq = 0.5, the NAMF has detection loss about 1.2 
dB w.r.t E-NAMF and 1.38 dB w.r.t GM-NAMF for the case of non-fluctuating 
target and about 2 dB w.r.t E-NAMF and 2.38 dB w.r.t GM-NAMF for the 
fluctuating model. This detection loss will increase for higher Pjj . For example, 
in the case of Swerling 0 target model at Pd = 0.9, NAMF has about 3 dB loss 
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Figure 4.9: FAP estimates. 



Figure 4.10: Gain estimates. 
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Figure 4.11: I g function estimates. 


GM-NAMF Detector, N = 64, L = 128, K = 50000 



Figure 4.12: Optimum biasing parameter estimates. 
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GM-NAMF Detector, N = 64, L = 128, K = 50000 



Figure 4.13: Threshold as function of FAP. 

GM-NAMF Detector, N = 64, L = 128, K = 50000 



Figure 4.14: Gain as function of FAP. 
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Swerling 0 target model, homogeneous Gaussian background, P fg = 10 6 



Figure 4.15: Pn in homogeneous case, Swerling 0 target. 


Swerling 1 target model, homogeneous Gaussian background, P =10 6 



Figure 4.16: Pd in homogeneous case, Swerling I target. 
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Figure 4.17: Pd in the presence of two interferers. Swerling 0 target. 


w.r.t E-NAMF for the same case and 5.5 dB loss w.r.t NAMF in homogeneous 
Gaussian background. It has also to be noticed that in the fluctuating target 
scenario the presence of nonhomogeneities produces a loss in term of Pd per¬ 
formance much more significant that in the nonfluctuating case, especially in 
the region of high SNR. From these results we conclude that in the presence 
of nonhomogeneities the E-NAMF and GM-NAMF detectors are much robust 
than the NAMF detector. Also, we compare the detection performance of the 
class of NAMF detectors with the AMF detector and its envelope and geometric 
mean variants. In Figure 4.19 and 4.20 these performance are shown respec¬ 
tively for the Swerling 0 and Swerling 1 target model. The nonhomogeneities 
considered here are the same as in the previous case. Properties of the two 
detector classes are summarized in Table 4.1. From these results, we conclude 
that the geometric mean and envelope detectors of the NAMF class have the 
best detection performance in the presence of interfering targets. 
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Swerling 1 target model, 2 interferers, P fa = 10 6 



Figure 4.18: Pd in the presence of two interferers. Swerling 1 target. 



Figure 4.19: Pd comparison in the presence of two interferers. Swerling 0 target. 
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Figure 4.20: Pp comparison in the presence of two interferers. Swerling 1 target. 


Table 4.1: Pp comparison of AMF and NAMF STAP detectors. 


Q 

AMF 

E-AMF 

GM-AMF 

NAMF 

E-NAMF 

GM-NAMF 

CFAR? 

Yes 

Yes 

Yes 

Yes 

Yes 

Yes 

homog. Pd? 

- 

e-loss 

e-loss 

- 

~ 0 loss 

~ 0 loss 

Robust? 

No 

Yes 

Yes 

No 

Yes 

Yes 
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4.6 Conclusion 


Envelope-law and geometric-mean variants of the NAMF detector have been 
introduced. Their respective thresholds have been determined for a range of 
FAPs. The detectors have been shown to have better performance in the pres¬ 
ence of outliers in the training data while maintaining almost equal performance 
as the standard square-law version in homogeneous Gaussian interference. To¬ 
gether with the variants proposed in Chapter 2, these detectors represent robust 
alternatives to conventional square-law processing. 
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Chapter 5 


Low-Rank STAP detectors 


5.1 Introduction 

In this chapter the STAP detector based on the low-rank approximation of the 
normalized adaptive matched filter (LRNAMF) is investigated for its perfor¬ 
mance. Being computationally efficient, the LRNAMF detector is a candidate 
for future implementation in STAP radar detection. Its FAP performance is 
analytically investigated here. 

In the second section the low-rank NMF detector is described. Subject to an 
approximation for the disturbance covariance matrix in a clutter dominated sce¬ 
nario, the FAP of the LRNMF detector is known via a simple formula, [16] and 
[17]. A brief description of the detector is reported here in order to summarize 
the low-rank approximation. The third section provides analytical derivation 
of the exact FAP of the LRNMF detector for data possessing an arbitrary co- 
variance matrix. The fourth section discusses a nominal covariance model for 
simulation and gives further analytical results. The fifth section discusses the 
LRNAMF detector and shows how the analytical results for the LRNMF detec¬ 
tor can be used to predict its performance. The sixth section provides simulation 
results followed by a conclusion. 

5.2 The LRNMF detector 

The low-rank approximation to the normalized matched filter is used for target 
detection in heterogeneous clutter scenarios. The background is assumed to 
consist of clutter plus white Gaussian noise. The binary hypothesis test can be 
written as 


Ho:x = d = c + n 

Hi:x = as + d = as + c + n (5-1) 
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where x is the primary data vector, c is the Gaussian clutter vector with co- 
variance matrix sR c with unknown level s and known structure, n denotes the 
additive white Gaussian noise vector with covariance matrix a 2 I, where I is the 
N x N identity matrix and the noise power a 2 is unknown, s is the steering 
vector and a is the unknown target amplitude. The vector d is used to represent 
the sum of the clutter and the white Gaussian noise. The covariance matrix of 
the disturbance d is 

R d = sR c + a 2 I (5.2) 

In many real cases the clutter covariance matrix R c has rank r which is less 
than N. This fact will be used to approximate the inverse of the disturbance 
covariance matrix which will be used in the NMF test. The matrix R^ can be 
expressed as 

R d = UDU t (5.3) 

where U is the matrix whose columns are the normalized eigenvectors of R^ and 
D is the diagonal matrix of the eigenvalues of R^. When R c has rank r <C N, 
then Rd can be rewritten as 


r N 

R d = ^(sXi + CT 2 )u ?; uj + ^ cr2 Wuj 

i—1 i=r +1 

The inverse covariance matrix is 

r N 

R.- 1 = J2(sXi + cr" 2 u iU J 

i—l i=r +1 

The previous expression can be rewritten as 

,t N 


= iz 




Z a 2u * u l 


i=r-\-1 


s\i 


1 f 1 0-2 t 

(72 ^i=l 11 a 2 (1 + 


i r 

?('-£ 

i= 1 


s\i 


a +m 


X Ui u- 

sXo \ 1 l 


(5.4) 


(5.5) 


(5.6) 


For a clutter-to-noise ratio (CNR) much greater than one, i.e. sA; a 2 , the 
inverse of the disturbance covariance matrix can be approximated as [18], 


where 




P = EL U ^ 


(5.7) 

(5.8) 
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is a rank r projection matrix formed with the r eigenvectors corresponding to 
the r dominant eigenvalues of R^. Using (5.7) in the NMF test of (3.1), we 
obtain the LRNMF test 


= 1st (I -P)x| 2 

LR - (st(I-P)s)(xt(I-P)x) 11 


(5.9) 


It can be noticed that the LRNMF test is invariant to the unknown clutter level 
s and to the noise power er 2 . Moreover, since (I — P) is also a projection matrix 
of rank (N — r), (I — P) 2 = (I — P), therefore we can define the transformed 
vectors si = (I — P)s, xi = (I — P)x and rewrite the test as 


Alr = 


l s l x i | 2 


H i 

i: 


(s{si)(xjxi) H, 


(5.10) 


It is important to observe that the low rank NMF test is the squared cosine 
of the angle between the transformed steering vector Si and the transformed 
data vector xi. Fast simulation can be performed using rotation of the primary 
data vector xi as described in Section 3.2.2 for the NMF detector; of course, 
since the two vectors in this case are not the same as for the NMF case, the IS 
performance will be different. 


5.2.1 FAP approximation: low clutter rank and high CNR 


The FAP for the LRNMF detector is derived in [17] and is given by 


a LRNMF 


(1 — r /)^-*— 1 


(5.11) 


A crucial point of the LRNMF detector is the clutter rank estimation. The 
threshold for the LRNMF test, for a fixed FAP, will depend on the rank r of the 
clutter covariance matrix, specifically it increases with increasing r. A technique 
for clutter rank estimation can be found in [19]. The case r = 0 coincides with 
the full rank NMF test, which is invariant to the white noise level. 


5.3 The LRNMF detector: arbitrary covariance 

R ( / 


In this section we derive two alternate expressions for the exact FAP of the 
LRNMF detector in (5.9). In particular, we do not assume that (I — P) in 
(5.7) whitens the primary vector X, as is required for the derivation of the FAP 
formula in (5.11). These exact forms are valid for any primary covariance matrix 
Rd. The first expression involves singular multivariate Gaussian distributions 
whereas the second does not. 

For ease of notation we denote the projection matrix (I — P) by Q, which is 
idempotent. The LRNMF test of (5.9) can then be written as 


Alr 


|stQ 2 x| 2 h .l 

(stQ 2 s)(xtQ 2 x) V 


(5.12) 


77 



Define the transformed vectors 


Then 


X! = Qx and Si = Qs 


Ri = -EjXxXj} = Q£{XX t }Q t = QR d Q f 


(5.13) 


(5.14) 


As Q is a singular matrix, the matrix Ri is also not of full rank. This fol¬ 
lows from the property of the rank of a product of matrices 1 . Hence Xj ~ 
5C7V w (0,Ri ), where SCAT indicates a singular multivariate Gaussian distribu¬ 
tion 2 . Using the above transformations the test becomes 

A = l s I x i | 2 = l s I x i ! 2 = l( s l/l|si||) xi | 2 = |ttxi | 2 

LR “ (S t 1 s 1 )(x t 1 x 1 ) || 8 i|| 2 ||x 1 ||2 || Xl p || Xl || 2 ^ j 

where t = si/||si|| is a unit vector. 


5.3.1 Exact FAP: using singular Gaussian distributions 

Define a unitary transformation G = Ht, such that the unit vector t-j has 
a single element equal to 1 and the remaining (N — 1) elements are zero, i.e., 
ti = (0,..., 0,1) T . The matrix H can be an Householder transformation matrix 


given by 

2uiP 

(5.18) 


H = I -^ 

l! u l! 

and 

u — t + | | e* 

V'i \ 

(5.19) 

where 

e, = [0,... ,0,1,0,... ,0]^ 

(5.20) 


1 rank(AB) < min( rank(A), rank(B)) 


2 A complex Gaussian vector x of TV components with covariance matrix R and mean vector 
fi is said to have a singular multivariate Gaussian distribution if rank(R) = p < TV. Using the 
factorization 


R = U 


D a 0 
0 0 


ut 


(5.15) 


the covariance matrix can be written as R = UiDaU), where U = [Ui,U 2 ] is an AT X TV 
orthogonal matrix, Ui is an TV X p column orthogonal matrix, and Da = diag(Ai . .... \ n ) with 
Ai > 0 for i = 1,... ,p. Define R~ = TJ ] D x 1 Tj). the generalized inverse of R = UiDaU 1 . 
The p.d.f. of x is given by 


/( x ) = * | exp{(x - /r) T R (x - n)} 

|Da| 


(5.16) 


where x lies on the p-dimensional linear subspace defined by 

L t (x-/i) = 0; L : TV X (TV -p), PR = 0, L f L = Ijv-p 
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has 1 as the ith element and zero elsewhere. The transformation by the matrix 
H will eliminate all the elements of t except the element f,;. Note that H is an 
Hermitian unitary matrix 

H -1 = = H (5.21) 

We will construct the H matrix in order to preserve only the IVth element of 
the vector t. The likelihood ratio in (5.17) now takes the form 


^LR 


|tlH Xl | 

xjxi 


(5.22) 


Now define 


y = H Xl 

R y = ^{yy 1 } = H^jxixljH 1 '= HRitf (5.23) 

where y ~ <SCA//v(0, R y ). The vector xi in the test can be replaced by 

x-i = HV = Hty = Hy (5.24) 


Then, the likelihood ratio of (5.22) and corresponding test can be rewritten as 

. |tiy| _ M 2 ^ 


^LR = 


jtjyl _ 

ytHtHy yty 


N 

E I Vi 

i=1 


^ V 


(5.25) 


where y,, i = 1 ,,N are the elements of the vector y. Now we can write the 
test in the form 


M 2 

Hi 

% 

Ho 

N-l 

^7(|s/jv| 2 + Y l^l 2 ) 

i= 1 


\y N \ 2 {l-rj) 

Hi 

i: 

Ho 

N-l 

v Y i^i 2 

i= 1 


M 2 

Hi 

Ho 

N-l 

Vo Y ^ 

i= 1 


M 

Hi 

Ho 


(5.26) 

where 770 = 77/1 — 77 . 

Define the (N — 1) vector y 

= (yi, ■ 

.. ,y N -i) T , Ml = \vn\, and 


Kv,y) = 

0 

M 

">l 2 ) 1/2 

(5.27) 
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Then, the FAP is given by 


e{p(u 1 > b( V , Y)|'y)} 

E{g{v, Y)} (5.28) 

and the g function can be expressed in integral form as 

/»00 

g(v, y) = / f(ui\y)dui (5.29) 

Jb(v,y) 

The conditional distribution of Y.y is also Gaussian and we have to determine 
the conditional mean and variance. The vector y can be written as y = [y; ijn] 
and the covariance matrix of y can be partitioned as 


“lrnmf 


A 


£{YYt} 

e{yy^} ‘ 


R-y ^yy N 1 

. E{Y n Y^} 

E{Y n Y^} _ 


-1 

* 
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(5.30) 


where Ry is the (N — 1) x (N — 1) singular covariance matrix of the vector y 

and er 2 is the variance of the random variable Ym. The conditional mean and 

y N 

variance of Yjy are 


N-l 

E{Y N \y} = Rj/ N yRy y = c T y = E c * ^ 

i =1 

var{F w |y} = a 2 yN - Ry N yRy Ryy„ (5.31) 

where c T = R Wv yRy 5 Cj, i = 1,..., N — 1 are the elements of the vector c, and 
R^ is the g-inverse of Ry. 

Conditioned on the vector Y, the random variable Yyr is Gaussian with 
mean and variance given in the expressions above. Then the random variable 
U\ = | Yn | has a Rice distribution with noncentrality parameter 


JV-l 


E JV —± 

. , CiVi 
1=1 


and parameter 


a 2 = var{Yjv|y}/2 


Hence the integral in (5.29) can be written as 

9 ( 1 7, y) = / ex P [ - (uj + s 2 )/2ct 2 ] Iq(u\S/ ex 2 ) dm 

Jb(v.v) ® 


(5.32) 

(5.33) 


(5.34) 
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and the FAP becomes 



(5.35) 


where /(y) is the density corresponding to the singular Gaussian distribution 
SCAf n-i(0, Ry). This is one expression for the FAP of the LR.NMF detector 
for arbitrary covariance matrix R<j of the primary data vector. 


5.3.2 Exact FAP: using nonsingular Gaussian distributions 

We start with the covariance matrix Ri of Xi ~ 5CA/"jv(0, Ri) having some 
rank say p < N in (5.14) and write its decomposition as 

Ri = QR d Q f = VAV* (5.36) 

where V is the matrix of eigenvectors corresponding to the eigenvalues {A,}^ 
in the diagonal matrix A, which can be written as 



(5.37) 


Define 


and 



(5.38) 


V p = upper left block of V G R pxp 
V N _ p = lower left block of V G R (JV “ p)xp 
Using these definitions we can write Ri as 

/ VpApVt VpApVLp \ 

V V^v-p ApVt Vjv-p Ap y\f_ p ) 


(5.39) 

(5.40) 


which also turns out to be the covariance matrix of a new vector Xi defined as 


1 Ap/ 2 W; W ~ CAf p (0,T) 

V N-P J F 
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and hence this new Xi is statistically identical to the vector Xi in the beginning 
of this subsection. This is a well known representation ([20]) for expressing a 
singular Gaussian vector in terms of a nonsingular Gaussian vector with inde¬ 
pendent components. 

Substituting this definition into (5.17) yields the likelihood ratio 

|tt Xl | 2 

Alr ll Y , 112 
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Ap /2 W 
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WtAp /2 

vt vL P ' 
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vp 
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N-p 


aV 2 w 


WtApW 


Now define 

Y = A* /2 W 

and the p x 1 unit vector 


(5.42) 

(5.43) 
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V"t 

p 


V 


N-p 


(5.44) 


Then Y ~ CJ\f p (Q. A p ) and the likelihood ratio becomes 


Alr 


|tjy| 2 

yty 


(5.45) 


Next we find a unitary transformation H (not to be confused with the unitary 
transformation H defined on page 78) such that 


ei = Hti = ( 1 • • 
Then the likelihood ratio can be written as 

|elHy | 2 


Alr = 


yty 


Defining 
leads to 

with Z ~ CA7 p (0, R s ) where 


Z = HY 


A _ l e I Z | 2 

Am = - 7 - 


z ' z 


FG = HApH 1 ’ 


(5.46) 

(5.47) 

(5.48) 

(5.49) 

(5.50) 
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The likelihood ratio test can therefore be written as 


N 2 ^ 

- P - ^ v 

E N 2 Ho 


where the z/s denote the components of z. Rearranging yields 


\ z i\ 
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Ho 
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2 ) 
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(5.51) 


(5.52) 


where ry 0 = 77/(1 — 77 ), which can be considered (somewhat loosely speaking) 
as the nonsingular version of the test in (5.26) on page 79. The rest of the 
derivation for FAP is similar to those in Section 5.3.1 but we go through the 
steps. 

Defining z = ( Z 2 • ■ ■ z p ) , the covariance matrix R 2 can be parti¬ 

tioned as 


R; 



zt 


£{|^i| 2 } E{Z 1 Zt} 
E{ZZ{} E{ ZZt} 


R 


R, 


Ri 


(5.53) 


The random variable Z\ is conditionally Gaussian with density f(z\ | z) cor¬ 
responding to the distribution CAf i(/i Cl a 2 ) where the conditional mean and 
variance /i c and a 2 are given by 

fj, c = E{Z 1 \z} = Rr^ R / 1 z 

= var (Z 1 | z) = of - K Zli R / 1 R^j (5.54) 

The random variable U = \Zi \ is then conditionally Rice with density 

f(u I z) = E exp T — (u 2 + s 2 )/2f7 2 l/o(us/cr 2 )| u>0 (5.55) 

CH 

where s = |/< c | and er 2 = er 2 /2. Defining 

1/2 

Kv, z) = (% l-il 2 ) ( 5 - 56 ) 

the FAP of the test in (5.52) can be expressed as 

a LRNMF = p{\Zi\ > ( ? ?0 5Z i=2 \Zi\ ) ) 

= e[p(u> 6 ( 77 , Z)|z)} 
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= E{g( V , Z)} (5.57) 

where 

/»oo 

< 7 ( 77 , z) = / exp [ — (u 2 + s 2 )/2(T 2 ]/ 0 ('us/cr 2 ) du (5.58) 

Jb(r), z) ^ 

using (5.55). Then, similarly to (5.35), the FAP can be written as 


“lrnmf J 
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r°° v 

1 exp r — (il 2 + s 2 )/2cr 2 l I 0 (us/a 2 ) /(z) du dz, 

6(77,z) ® 

( 5 . 59 ) 


where /(z) corresponds to the (nonsingular) density of CA/’ p _i(0, Rz). This is 
the second expression for the exact FAP of the LRNMF detector for arbitrary 
covariance matrix Rd of the primary data vector. 


5.4 Nominal statistical model for simulation and 
threshold setting 


A nominal statistical model for the radar returns in the (target-free) primary 
data vector is required in order to be able to specify a covariance matrix Rd 
under the Gaussian assumption. Such a model can then be used to derive thresh¬ 
old settings for desired FAP values for the detector either through simulation 
or through direct computation of FAP expressions. 

For the disturbance covariance Rd = sR c + cr 2 I in (5.2) we start again with 
the low clutter-rank approximation of (5.4) on page 76 as 

r N 

R d = ^(sA j; + cr 2 ) Ujiij + ^2 

i— 1 i=r +1 

= UDU t (5.60) 

and identify the matrix of eigenvalues D as 

D = diag ( sAi + a 2 , ..., s\ r + a 2 , a 2 ,..., cr 2 ) 

= cr 2 diag ( 1 + sXi/a 2 ,... , 1 + sA r /er 2 ,1,..., 1 ) (5.61) 

N-r 


Now, choice of the unitary matrix U of eigenvectors will determine the covari¬ 
ance matrix R^. The simplest choice that can be made (at the risk of seeming 
somewhat subjective) is U = I. In such a case, Rd = D. That is 


( 1 + sAi/tr 2 0 

0 


Rd = D = 


1 + s\ r /a 2 


0 \ 


(5.62) 


. 01 / 
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It turns out that this choice is not without some practical significance. In a 
special case, using the high clutter-to-noise ratio assumption 


s 4 »i 

a z 


with 


a 2 = 1; A, = 1, i = 1,... ,r 


leads to a covariance matrix of the form 


Rd — D — 


( s 0 
0 s 


Vo 


' ° ^ 

1 0 

0 l) 


(5.63) 


A plot of the eigenspectrum of R^ estimated from Gaussian data generated 
according to the above R^ is shown in Figure 5.1. This eigenspectrum shows 
a close resemblance to the eigenspectrum estimated from KASSPER data, [17]. 
We note of course that this property of our model is due to the shape of the 
matrix of eigenvalues D (rather than the particular choice U = I). 

A consequence of the structure of the nominal covariance matrix R^ in (5.62) 
is described in the following subsection. 


5.4.1 Exact FAP of LRNMF detector: nominal 

We calculate here the FAP of the LRNMF detector when the covariance matrix 
Rd of the primary data vector is as in (5.62). The assumption of high clutter- 
to-noise ratio is not used. Therefore, the only restriction is on the shape of R<j. 
As Rd is diagonal, it follows that its eigenvector matrix U are composed of the 
columns of I and hence, for a given r, the projection matrix Q is given by 


r 

Q = I - Y UiU i 

i— 1 

N-r 

= Y U i U i 

( 0 . . . . 0 \ 


Vo 



. 0 
0 1 j 


(5.64) 


which is of rank N — r. Applying this and R^ in (5.62) to (5.36) yields 

R x = VAV t 
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(5.65) 


This implies, from (5.37), that 


ll 0 
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0 \ 

0 
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(5.66) 


p = N — r, and thus A p in (5.38) is A p = lAr_ r . Using this in (5.50) yields 
R 2 = Uv— r . Therefore the test in (5.52) represents a simple CA-CFAR test 
with FAP given by 


ft LRNMF 


(i - v ) N - r - 1 


which coincides with the expression for the LRNMF detector in (5.11) on page 77. 

Therefore we can conclude that whereas the above formula is an approxima¬ 
tion (albeit a good one) for the FAP under the assumptions of low clutter-rank 
and high clutter-to-noise ratios, it is an exactitude if the covariance matrix of 
the data possesses the structure in (5.62). In the latter case the clutter-to-noise 
ratio does not matter. 


5.5 The LRNAMF detector 


The low-rank normalized adaptive filter (LRNAMF) detector is described by 


A 


LR-A = 




(st(I-P)s)(xt(I-P)x) 


H i 

H 0 


(5.67) 


where 

p = ^2 , uiu i ( 5 - 68 ) 

and {uj}^ are eigenvectors of the covariance matrix estimate obtained from 
the secondary vectors with r being the estimated rank of the clutter component 
of the disturbance. 


5.5.1 FAP approximation: low clutter rank and high CNR 

The main difference between the LRNMF and LRNAMF detectors (the former 
given in (5.9)), for a fixed value of r used in both detectors, lies in the set of 
eigenvectors {u,}( v used to compose the matrix P (or P). However, for the 
low-rank clutter model of (5.4) and subject to the high clutter-to-noise ratio 
approximation of (5.7), the FAP performance of the LRNMF detector is given 
by the expression of (5.11) and this is independent of the eigenvectors used in the 
detector. This implies that the LRNMF detector is CFAR with respect to the 
structure of covariance matrices R^ of data that arise from models satisfying 
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these assumptions. Therefore, for data having the same eigenspectrum as in 
the above model, and for the same value of r used in both detectors, the FAP 
performance of the LRNAMF detector will be, largely, equal to that of the 
LRNMF detector and given by the same formula for each realization of sample 
covariance matrix R^. The LRNAMF detector will also be approximately CFAR 
under the above assumptions on clutter rank and power. Nevertheless, it must 
be emphasized that the FAP expression of (5.11) is an approximation. 

The above conclusion is borne out by simulation results. 

5.5.2 Exact FAP: arbitrary covariance 

If the model assumptions do not hold, specifically in cases where the CNR is 
not high, then the expressions for exact FAP of the LRNMF detector derived 
in (5.35) and (5.59) can be used to determine the exact FAP of the LRNAMF 
detector. Whereas the influence of the secondary vectors on the LRNAMF 
detector is felt through the eigenvectors in P, via singular value decomposition 
of the covariance matrix estimate R,/, the FAP is affected by the covariance 
matrix R z in the expression (5.59). This matrix is derived from the eigenvectors 
in P through a series of transformation described in Sections 5.3 and 5.3.2 and 
is, strictly speaking, a random matrix whose statistical properties are related to 
those of the secondary vectors. Therefore the FAP expression for the LRNMF 
detector in (5.59) can be considered as a conditional probability expression; the 
average of this over the statistics of the covariance matrix R z will produce the 
exact FAP performance of the LRNAMF detector. Although it does not appear 
to be analytically tractable, this fact can be formally expressed by 

Q LR.NAMF = -®'R*{ Q! 3(Rz)} 

= £r z {£{< 7(77,Z)|R z }} (5.69) 


with g(r], z) given in (5.58). 


5.6 Simulations for LR detectors 

We describe here the simulation results obtained for FAP and threshold estima¬ 
tion for both LRNMF and LRNAMF detectors for various assumptions on the 
data covariance matrix R^. The two FAP expressions, (5.35) and (5.59), have 
been used to set up estimators for FAP. The estimator for the LRNMF detector 
can be written as 

1 K 

“lrnmf = Z); Z ~ C7V p _i(0, R z ) (5.70) 

1 
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For simplicity we have not used any IS. The FAP estimator for the LRNAMF 
detector is configured as 


“lrnamf 


1 K 2 . K 1 ... 

* i 1 i 

Z ~ C7V p _i(0, Rz), Rz ~ /(secondary vectors) (5.71) 


where the inner simulation is conditioned on the covariance matrix R^, which 
depends on the data covariance R^ in place. 

Detector thresholds have been estimated using the inverse method. 


5.6.1 models used 

Three different data covariance models have been used to study the perfor¬ 
mances of LRNMF and LRNAMF detectors for FAP as well as detection prob¬ 
ability. The models are 

Case 1 : In this model, the disturbance is assumed to be white Gaussian noise 
which ideally has a flat eigenspectrum. 

Case 2 In this model we use an eigenspectrum with eigenvalues ranging from 
10 4 to 0.1, covering the same range as Model 1, but linearly (in dB) with 
eigenvalue number. 

Case 3 This is the nominal disturbance model which has the matrix of eigen¬ 
values described in (5.63), and corresponds to the case of low clutter-rank 
and high CNR. 

Data generated using these models have estimated eigenspectra shown in Fig¬ 
ure 5.1. Wherever detection probability results are provided for a certain model, 
the detector thresholds used have been estimated to provide the stated FAPs 
for that model. 


5.6.2 Results for FAP 

The first result for FAP of the LRNMF detector is shown in Figure 5.2 for 
the eigenspectrum of Case 3. As expected, the simulation results using the 
//-method coincide with the formula of (5.11), indicated by stars in the figure. 
Corresponding simulation gains for the (/-method are shown in Figure 5.3. 

The next figure, Figure 5.4 shows results for FAP versus thresholds of the 
LRNMF and LRNAMF detectors for the three eigenspectra models. Lack of 
CFAR is clearly observed for both the LRNMF and LRNAMF detectors with 
change in the covariance structure of the data. However, a surprising obser¬ 
vation is that for the same data eigenspectrum, the LRNMF and LRNAMF 



Eigenspectrum 



Figure 5.1: Estimated eigenspectrum for the 3 data models. Clutter rank is 
r = 33. 


detectors match very closely in terms of achieved FAP values for given thresh¬ 
olds 3 . A plausible explanation for this behaviour lies in the following arguments. 
Firstly, LR detectors seem to ignore some of the information available in the 
secondary vectors, depending solely on the eigenvectors obtained through spec¬ 
tral decomposition and not on the actual estimated eigenvalues. Furthermore, 
for a sufficiently large number L of secondary vectors, the estimated eigenspec¬ 
trum resembles the actual eigenspectrum in shape. This fact may manifest 
itself through the eigenvectors of the spectral decomposition. We emphasize 
that although the match between LRNMF and LRNAMF detectors appears 
close graphically, numerical estimates obtained for the two detectors do differ. 

Estimated simulation gains for the LRNAMF detector are shown in Fig¬ 
ure 5.5. 

5.6.3 Results for detection probability 

Various detection probability results are shown in Figures 5.6 - 5.9. Figure 5.6 
shows that the detection probability of the LRNAMF detector is almost un¬ 
changed from the homogeneous case in the presence of interfering targets in the 

'’The exception appears to be Case 2. However this is not a departure from the above 
observation as the number of simulation trials for FAPs below 2 X 10 —4 were insufficient, 
according to gain estimates that were made. Hence more simulations are required to make a 
definitive conclusion for this eigenspectrum. 
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Figure 5.2: Comparison of the ^-method estimator and formula (5.11) for FAP. 
N = 128, r = 33. Eigenspectrum is for Case 3. 


secondary data, for both Swerling 0 and 1 target models. Comparisons with 
other detectors are shown in Figure 5.7 and Figure 5.8 for the eigenspectrum of 
Case 3. Comparing the performance of the adaptive detectors (LRNAMF and 
NAMF) for this model, the LR version performs better than the NAMF detector 
for both Swerling 0 and Swerling 1 target models, not only in the homogeneous 
case but also in the presence of interfering targets contaminating the secondary 
vectors. On the other hand when the background is as in Case 1 (i.e., clutter has 
rank zero), the NAMF detector performs better than the LRNAMF detector in 
the presence of interfering targets, as shown in Figure 5.9 and Figure 5.10. In 
homogeneous background the NAMF detector performs better for the Swerling 
0 target model and the two detectors perform the same for Swerling 1 target 
model. 
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Figure 5.3: Gain of the ^-method for N = 128, r = 33. Eigenspectrum is for 
Case 3. 



LR Detectors, N = 64, L = 128, r = 33 


case 2 LRNMF 
case 3 LRNMF 
case 1 LRNMF 
case 1 LRNAMF 
case 2 LRNAMF 
case 3 LRNAMF 


10“ 91 — 
0.05 


0.25 


0.35 


Figure 5.4: Thresholds versus FAP for three eigenspectrum models. 
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LRNAMF, N = 64, L = 128, r = 33 



Figure 5.5: Gain of the ^-method for N = 64, L = 128, r = 33. Eigenspectrum 
is for Case 3. 


LRNAMF detector, N=64,L = 128, r = 33, P, =1 O' 6 

fa 



Figure 5.6: Detection probability for LRNAMF detector, N = 64, L — 128. 
Eigenspectrum is for Case 3. 
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Swerlin 0, N=64,L = 128, r = 33, P, =1 O’ 6 

fa 



Figure 5.7: Detection probability comparisons, N = 64, L = 128. Eigenspec- 
trum is for Case 3. 


Swerling 1, N=64,L= 128, r = 33, P =1 O' 6 



Figure 5.8: Detection probability comparisons, N = 64, L = 128. Eigenspec- 
trum is for Case 3. 
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Swerling 0, R d = I, P fa =10 -6 



Figure 5.9: Detection probability comparisons (Swerling 0), N = 64, L = 128. 
Eigenspectrum is for Case 1. 


Swerling 1, Ft = I, P fa =10" 6 



Figure 5.10: Detection probability comparisons (Swerling 1), N = 64, L = 128. 
Eigenspectrum is for Case 1. 
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5.7 Conclusion 


In this chapter some theoretical investigations have been made into the perfor¬ 
mance of low-rank STAP detectors. Principally, we have been able to character¬ 
ize the FAP performance of LRNAMF detectors in terms of detection thresholds 
and disturbance backgrounds. Several detection probability simulation results 
have also been obtained. Some of these results are surprising and might seem 
somewhat counterintuitive. We have not been able to find complete justifica¬ 
tions for these observations and clearly further work is required. 
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Importance sampling for characterizing STAP 

detectors 

Raj an Srinivasan and Muralidhar Rangaswamy 


Abstract —This paper describes the development of adaptive 
importance sampling techniques for estimating false alarm prob¬ 
abilities of detectors that use space-time adaptive processing 
(STAP) algorithms. Fast simulation using importance sampling 
methods has been notably successful in the study of conventional 
constant false alarm rate (CFAR) radar detectors, and in several 
other applications. The principal objectives here are to examine 
the viability of using these methods for STAP detectors, develop 
them into powerful analysis and design algorithms and, in the 
long term, use them for synthesizing novel detection structures. 
The adaptive matched filter (AMF) detector has been analyzed 
successfully using fast simulation. Of two biasing methods con¬ 
sidered, one is implemented and shown to yield good results. 
The important problem of detector threshold determination is 
also addressed, with matching outcome. As an illustration of 
the power of these methods, two variants of the square-law 
AMF detector that are thought to be robust under heterogeneous 
clutter conditions have also been successfully investigated. These 
are the envelope-law and geometric-mean STAP detectors. Their 
CFAR property is established and performance evaluated. It 
turns out the variants have detection performances better than 
those of the AMF detector for training data contaminated by 
interferers. In summary, the work reported here paves the way 
for development of advanced estimation techniques that can 
facilitate design of powerful and robust detection algorithms. 

Keywords STAP, AMF, Importance Sampling (IS), CFAR, 
GM-STAP, E-AMF, Homogeneous clutter. Complex Gaussian 
PDF 

I. Introduction 

Estimation of false alarm probabilities of detection algo¬ 
rithms that employ space-time processing is examined here 
using forced Monte Carlo (MC) or importance sampling (IS) 
simulation. Space-time adaptive processing algorithms are 
of much importance for radar detection, [1], [2]. They are 
notoriously intensive from a computational point of view, with 
the more advanced (and robust) ones being also analytically 
difficult to quantify. Therefore, it is appropriate to attempt to 
develop fast simulation methods that could be used in their 
analysis and design. 

The work of the first author was supported by the European Office of 
Aerospace Research and Development, London, UK, in collaboration with 
the Air Force Office of Scientific Research under Award No. FA8655-04-1- 
3025. Portions of this paper were presented at the IEEE International Radar 
Conference, Washington, DC, USA, May 2005. 

R. Srinivasan is with the Telecommunication Engineering Group, Uni¬ 
versity of Twente, PO 217, 7500 AE Enschede, The Netherlands (email: 
r. srinivasan @ ewi. utwente. nl). 

M. Rangaswamy is with the Air Force Research Lab¬ 
oratory Sensors Directorate, Hanscom AFB, MA, USA 
(email: Muralidhar. Rangaswamy @ hanscom. af. mil). 


In the following we use lessons learnt from developing 
IS techniques for characterizing conventional CFAR detec¬ 
tors and describe an experiment in applying them to STAP 
detection. The starting point of this effort is the celebrated 
AMF detector derived in [2] and which represents the array 
version of the workhorse cell averaging (CA) CFAR detector 
for conventional radar signal processing algorithms. The false 
alarm probability (FAP) performance of the AMF detector is 
known in integral form under certain conditions and can be 
numerically computed to any desired accuracy. Thus it forms a 
suitable basis for validating our simulation experiments. Two 
specific IS methods (described in the sequel) are presented 
and the better (and also easier) one is implemented. As a 
demonstration of the power of IS methods, the envelope- 
law and geometric mean detectors are presented and analyzed 
using fast simulation. These are variants of the AMF detector 
and their FAPs are not known in analytical form. Their CFAR 
property is established and FAP behaviour characterized using 
IS. A brief comparison of detection performances of all 3 
detectors is made. In the long term, our aim is to develop 
and apply these fast simulation techniques to modern STAP 
detection algorithms (see [3] - [5], and references therein). An 
important goal in this context is to devise IS biasing methods 
that result in simulation times that grow slowly with decreasing 
FAPs. This remains an open problem. 

As well known now, IS is the chief simulation methodology 
for rare-event estimation. It is an enduring method that has 
had success in several areas of science and engineering, [6]. 
Briefly, IS works by biasing original probability distributions 
in ways that accelerate the occurrences of rare events, con¬ 
ducting simulations with these new distributions, and then 
compensating the obtained results for the changes made. The 
principal consequence is that unbiased probability estimates 
with low variances are obtained quickly. The main task in IS 
therefore is determination of good simulation distributions for 
an application. Simulations performed using such distributions 
can provide enormous speed-ups and, if applied successfully, 
simulation lengths needed to estimate very low probabilities 
can become (only) weakly dependent on the actual probabili¬ 
ties. When these probabilities satisfy a large deviations princi¬ 
ple, then several asymptotic results are available for devising 
simulation distributions, [6]. For most detection applications 
however, it appears that adaptive methods which attempt to 
minimize error variances might be better suited. There have 
been some recent attempts in the literature, for example [7], 
[8], to apply IS for FAP estimation of CFAR detectors with 
varying degrees of success. Our work uses results developed 
in [9] and [10], with an adaptive implementation of the so 
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called g-method. Recent theoretical and application work on 
this method can be found in [11] and [12], The groundwork 
for the present results was laid down in [15] which, to our 
knowledge, is the first article on the topic of using IS for 
studying STAP detection algorithms. 

During the conduct of simulations reported herein, some 
issues concerning the adaptive IS algorithms used arise, and 
these are discussed briefly. The positive outcome of the 
methods used is that excellent match with numerical results is 
obtained. The succeeding sections provide a brief statement 
of the detector tests, how IS biasing can be performed to 
hasten false alarm events, description of the g-method which 
is a conditional IS technique developed originally for studying 
sums of random variables [9], how inverse IS can be used to 
estimate (and choose) detector thresholds, a methodology for 
using the fast adaptive algorithms, simulation results, and a 
conclusion. 


which can be used to numerically determine the threshold 
setting for a desired FAR As shown in [2], the test in (1) 
can be rewritten as 

|s^R _1 x| 2 ^ ?7 S^R -1 s 

Ho 

= (3) 
n i=1 

This is in the form of a vector (or, array) version of the 
usual CA-CFAR test. The LHS is a square law detector, 
being the output of a matched filter (matched to the direction 
s in which the array is steered) for incoherent detection 
using the so-called sample matrix inversion (SMI) beamformer 
weights If s. The RHS represents a cell averaging term. 
Further details on these issues can be found in the references 
mentioned above. 


II. STAP DETECTOR VARIANTS 

In a radar system consisting of a linear array of N s antenna 
elements a burst of N t pulses is transmitted, resulting in a 
return of N s N t = N samples in each range gate. The N 
complex samples are arranged in an N x 1 column vector 
and there are L + 1 such returns. The range gate return to 
be tested for the presence of target is called the primary data 
and is denoted by x while the remaining vectors are known 
as the training (or secondary) data, assumed to be free of 
target signal, and denoted by x(7), l = 1 A target 

return is modelled as consisting of an N x 1 steering vector s 
with an unknown complex amplitude in addition to clutter, 
interference, and noise. The primary and secondary data 
vectors are assumed to be jointly independent and zero-mean, 
spherically symmetric complex Gaussian, sharing the N x N 
covariance matrix R = /i{XX ; }, where the superscript f 
denotes complex conjugate transpose. 


A. The AMF (square-law) STAP detector 

Under these assumptions the AMF detection test, as ob¬ 
tained in [2], is given by 


IstR-ixl 2 h , 

—- € 7 1 

stR _1 s H 0 


(1) 


where 

R= j-^x(l)x(l) f 

n ;=i 

is the estimated covariance matrix of x based on the secondary 
data (also referred to as sample matrix), and rj is a threshold 
used to set the FAP at some desired level. This test has the 
CFAR property. The FAP a of the test is known to be given 
by 


LI 


1 X L-N +1 


(1 - x) 


N—2 


(L — N + 1)!((V — 2)! J 0 ( l + r,x/L ) L ~ N + 1 


dx 


( 2 ) 


B. The AMF (envelope-law) STAP detector 

In contrast to square-law detectors, it is known to radar 
engineers that envelope detectors possess some robustness 
properties in terms of detection performance when the training 
data is contaminated by outliers or inhomogeneities. Accord¬ 
ingly, the envelope-law STAP version of the AMF (abbreviated 
hereafter as E-AMF) detector is proposed here and its detec¬ 
tion properties evaluated. It takes the form 

H L 

| S tR-i X | i dK V'|stR- 1 x(0| (4) 

H ° L t! 

where r] E denotes the threshold multiplier. An analytical 
expression or approximation for the FAP of this detector is 
not known. 


C. The geometric-mean STAP detector 

Also proposed here is the geometric-mean STAP version of 
the AMF detector (abbreviated as GM-STAP). Conventional 
CFAR detectors that calculate the geometric means of the 
range cells in the CFAR window (usually referred to as log-t 
detectors) are known to have robustness properties. The STAP 
variant 1 takes the form 


sR 


-i. 


H i 
Ho 


1/L 




ix (oi 


(5) 


Kl =1 


where ?/ r; denotes the threshold multiplier. Note that the 
GM-based square- and envelope-law versions are identical 
(except for a trivial squaring of threshold multiplier). For 
the sake of completeness, the value of the multiplier as 
the number of training vectors L —> oo is calculated. In 
such a case, R 1 R 1 , stR X s^R *X, 

and R X(7) s^R -1 X(7) in the absence of tar¬ 

get. As sTR _1 X and sTR _1 X(7) are i.i.d. and distributed 
as CA/”i(0,s^R^ 1 s), it follows that the squared envelopes 
| s^R _1 X| 2 and |s^R _1 X(/)| 2 are exponentially distributed. 


lr This test was actually suggested in [16] but its performance was not 
evaluated. 
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each with mean s'R hs. Using convergence arguments based 
on continuity, [14], it is straightforward to show that the 
asymptotic FAP of the GM detector is exp (—r / 2 exp (— 7 )) 
where 7 is the Euler constant. Hence the (asymptotic) thresh¬ 
old rj G required to provide a FAP of 1CU 6 is 4.9605. A similar 
calculation can be easily made for the E-AMF detector. 

Both the E-AMF and GM-STAP detectors have the CFAR 
property under the assumption of homogeneous Gaussian 
characterization as stated in the beginning of this section; their 
FAPs are independent of the true target-free covariance R. 
This is shown in Appendix I. 

III. FAP ESTIMATION USING IS 

This section describes procedures to quickly estimate the 
FAPs and threshold multipliers of STAP detectors using IS. 
The standard (square-law) AMF detector is used for ex¬ 
position. The detector variants are handled using parallel 
arguments. 

The actual form of biasing considered here is simple scaling. 
While other techniques (superior in terms of simulation gains) 
can be devised, scaling is easy to implement and yields 
conservative but reasonable results. 

A. Square-law AMF detection 

Two strategies to estimate FAPs using IS are two- 
dimensional ( 2 -d) biasing and the conditional (/-method pro¬ 
cedure, described in this section. The 2-d biasing scheme 
parallels the approach used in previous works on conventional 
CFAR algorithms, as cited in the introduction. It can be useful 
in situations wherein the more powerful (-/-method might be 
difficult to apply. 

1) 2-d biasing: To estimate FAP using IS, we make the 
following observations. Suppose each complex sample of a 
secondary vector is scaled by a real (and positive) number 
f? 1 / 2 . This has the effect of scaling the covariance matrix 
estimate FI by 6. Therefore, as far as the covariance estimate 
is concerned, both sides of the test in (3) remain unaffected 
by the scaling. However, each secondary vector being scaled 
by 0 1 / 2 results in an additional scaling of the RHS by 9. 
Hence choosing 9 less than unity will have the effect of 
compressing the density function of the random threshold of 
the test. Further, a scaling of each complex component of the 
primary vector by a real positive a 1 / 2 will achieve a scaling 
of the LHS of the test by a. Thus, choosing a larger and 9 
smaller than unity will achieve an increase in the frequency 
of occurrence of false alarm events during simulation. The IS 
optimization problem will be a two-parameter one. 

Denoting by A the false alarm event in (3), the unbiased IS 
estimator in an i.i.d. (independent and identically distributed) 
simulation can be expressed as 

1 K 

3 = -^[lM)W(X,X L ;a,0)]W; X,X L ~/* (6) 

i=1 

where K is the length (or number of trials) of the IS simu¬ 
lation, l(-) denotes the indicator, X^ = (X(l),... ,X(L)) r , 
and the notation ~ /* means that all random variables are 
drawn from biased distributions. The weighting function W 


is described below. In setting up the joint densities of X and 
Xi, we use the fact that the FAP of the AMF has the CFAR 
property and is independent of the true covariance matrix R. 
This is so under the assumption of Gaussian distributions for 
the data. In such a case, the simulation of the AMF test can be 
carried out for data possessing a diagonal covariance matrix 
I, denoting the N x N identity matrix. Therefore, primary 
and secondary data can be generated as complex vectors with 
independent components. The unbiased joint densities are 

cW -V I x (!)t x (() 


/(x) = 


tN 


so that 


/(x,x L )= 


and /(x L ) = — 

g-x+x-Ef x(i) f x(0 


LN 


n (L+l)N 

With scaling performed as described above, the biased joint 
density takes the form 

rt x-s Ei x(I) + x(t) 


— x x- 


/*(x,x L ) = 


n (L+l)N a NQLN 

resulting in the weighting function 

/(x,x t ) 


W (X, Xij a, 9) = 


/*(x,x L ) 

= Ca N 9 LN e A ' a e B ' s 


(7) 


where 


A = x^x, B = x(Z)^x(Z), and C = 


= P ~{A+B) 


( 8 ) 


The variance of the IS estimator a can be expressed as 

var a= ^[I(u) - a 2 } (9) 

with v denoting the vector biasing parameter ( a,9) T £ 
[l,oo) x (0,1]. The quantity / is given by 

I{y) = E i ,{l(A)W 2 (X,X L ;u)} 

= E{l(A)W(X,X L ;u)} (10) 


where the expectation E+ proceeds over biased distributions. 
Minimization of vara with respect to the biasing parameters 
is equivalent to minimization of I and is described in Ap¬ 
pendix II. 

2) The g-method estimator: This method exploits knowl¬ 
edge of underlying (input) distributions more effectively, yield¬ 
ing a more powerful estimator. Additional advantages are that 
only a scalar parameter optimization problem needs to be 
tackled and the inverse IS problem (for threshold optimization 
or selection) can be easily solved. The FAP can be written as 

a = P(A) 

= £’{P(|s t R _ 1 X | 2 > r/s t R~‘ 1 s | X L ,H 0 )} 

= E{g(X L )} (11) 

where (/(Xj) denotes the conditional probability in the second 
step. The conditioning implies that the covariance matrix 
estimate R is given. We proceed to estimate a using the form 
in the third step above. 
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With the conditioning in mind it is easy to show, assuming 
that X is rotationally invariant and Gaussian, that the random 
variable slR _1 X = w^X is distributed as CJ\f i(0, w*Rw) 
with independent real and imaginary components, and the 
weight vector w = 11 '.s. Hence the random variable Y = 
|s^R -1 X | 2 is exponential and has the density function 

-j/w t Rw 

f(y\X L ,H 0 ) = — --, y> 0 

wTRw 

Therefore 

g(X L ) = P(Y > tys^R _1 s | X L ,H 0 ) 


where 


stR-is 

wfRw 


Note that if R = R, then r/fX/J = e~ v and this is the FAP of 
the AMF when the covariance matrix is known. As discussed 
before, we run simulations with homogeneous data possessing 
an identity covariance matrix, that is, with R = I. The g- 
method IS estimator then takes the form 

= ^J2[g(X L )W(X L -e)}^ 

i =1 

1 K 

= ( 12 ) 

i =1 


with D given by 


D 


slR-^ 

S lR _1 s 

st(R-i)2 s 


(13) 


Now, choosing the (single) biasing parameter 8 < 1 pro¬ 
duces a decrease in D , thereby causing a higher frequency of 
occurrence of the false alarm event or, more appropriately in 
this case, a larger value of the (/-function. Note that use of 
the (/-method obviates the need to bias primary data vectors. 
Determination of a good value of 8 proceeds as described in 
Appendix II. The weighting function is simply 


W(x L ; 8) = QLN e -(l-l/9)B ( 14 ) 


which can be deduced from (7) by setting a = 1. The variance 
of this estimator is given by 

vara g = ^p[I g {0) - a 2 ] (15) 

A 

where 


I g (8) = E*{g 2 {-K L )W 2 {X L -8)} 

= E{g 2 (X L )W(X L -8)} (16) 


The minimization of (an estimate of) I g with respect to the 
biasing parameter 8 is carried out using the recursion 


@m +1 — ( 1 /J 


r fgifim) 


which is just a one-dimensional version of (35) in Appendix II. 
The reader is again referred to this appendix for definitions of 
quantities in the summands of the following estimators for I g 
and its derivatives (with respect to 8) that are used in (17). 
These are respectively given by 

im = ^E^ 2 ( x ^ 2 ( x ^)] w 

i=1 

? g (0) = ^Y J [g 2 {X L )W{^8)W e {X L -d)}^ 

i= 1 

1 K 

= -^[g 2 (X L )W(X L ;0)lH ee (X L ;0)]« 

i=l 

(18) 

where X^ ~ /* in all 3 estimators above. 

3) Threshold determination: The converse problem, namely 
that of finding by fast simulation the value of detector thresh¬ 
old 77 satisfying a prescribed FAP, is an important task that can 
be readily carried out using the estimator of (12). The genesis 
of the method lies in the so called inverse IS problem first 
formulated and solved in [9]. The solution is to minimize a 
“squared performance error” stochastic objective function 

J(v) = [a g (v) - “of 

where a 0 is a desired FAP. A typical example is shown in 
Figure 1 for an E-AMF detector. All detection algorithms that 
involve a threshold crossing will possess objective functions 
that have the general behaviour shown for a 0 < 1 , assuming 
of course that the FAP estimate is (in a stochastic sense) 
monotone in the threshold //. Minimization of J with respect 
to 77 is carried out by the algorithm 

_ . c a o ~ OtgiVm) _ 1 0 no\ 

T]m+1 — Vm T Or) . — - 5 TO — 1,2,... (19) 

Ot'giVm) 

where 5 V is a step-size parameter and the derivative estimator 
in the denominator is given by 

„ , k 

a' g (n m ) = - 1{ J2^ De ~ vDw ( x P e ^ i) ^ ( 20 ) 

i =1 

with the prime indicating derivative with respect to 77 . Note 
that this derivative estimator actually estimates (negative of) 
the probability density function of the AMF statistic on the 
left hand side of ( 1 ) under the hypothesis that target is absent. 
The threshold-finding algorithm of (19) is a key component 
of the fast simulation methodology. 

4) Simulation gain: A useful measure of the effectiveness 
of any IS procedure is the simulation gain T. It is the ratio 
of simulation lengths required by conventional MC and IS 
estimators to achieve the same error variance. The variance of 
an MC estimator (W = 1) is given by (a — a 2 )/K M c for a 
simulation length Kmc. Hence, using (9), the IS gain for the 
estimator of ( 6 ) is 

P _ Kmc _ ol — a 2 
~ ~~K~ ~ I(u) - a 2 


TO — 1,2,... 


(17) 


( 21 ) 
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with I(y) defined in (10). Similarly, using (15), the (/-method 
estimator of (12) has gain l' r/ over MC estimation given by 


r 


a 


a — a 2 
I g (0) - a 2 


( 22 ) 


where I g (9) is given by (16). Note that setting 9 = 1 in the 
above provides gain of the (/-method estimator without IS. 
The estimator always has a smaller variance and consequently, 
T g > T. Comparing (21) and (22), it is sufficient to show that 
I g < I to prove this. This is done in Appendix III. 

As seen in the sequel, estimation of gain during simulation 
plays an important role in mechanizing an IS methodology 
that allows rare events to be studied quickly. In particular, T g 
is estimated as 


r s = 


a Q — 


m-- 


(23) 


wherein the required estimates are given in (12) and (18). The 
denominator of the above equation is related to an estimate of 
the IS variance, denoted as var a g . Substituting (12) and (18) 
into (23) and performing a little algebra yields 


var a g 



»= 1 3=1 

> 0 (with probability 1) (24) 

where y t = [</(X L ) W(X L ; 0)]». The variance estimate is 
asymptotically unbiased and it is easy to show, using (15) and 
(24), that var a g —vara g —> 0 as K —■> oo, for any fixed 9. 

If the IS estimator of (6) is used instead of (12), then 
j/,; = [l(-4) W(X,Xl\ a, 0)]W and the inequality in (24) will 
not be strict, leading to a non-zero probability of instability 
in gain estimation. One way of trying to prevent this from 
happening is to overbias the simulation, but this can result in 
underestimation of the false alarm probability. Such instabili¬ 
ties do not occur with the (/-method estimator and this is one 
of its implementation advantages. 


B. E-AMF and GM-STAP detection 

False alarm probability IS estimators for these two detector 
variants are set up in completely parallel fashion, starting with 
the expression in (12). The only difference is in the definition 
of the quantity D (given in (13) for the AMF detector), 
while the threshold multipliers get squared. This follows by 
observing from (3), (4), and (5) that all 3 detectors studied 
here differ, primarily, on the right hand sides of their respective 
tests. Consequently, denoting by De and Dq respectively the 
D values for the two variants, this leads to 

_ (rEf=.i |s t R- 1 x(Z )|) 2 

De = --—^- 

st (R _1 ) 2 s 

and 

_ (nf = iis t R- i x(oi 2 ) 1/i 

Dq = -^- 

st(R-i)2 s 


The corresponding FAP estimators are obtained by us¬ 
ing the above definitions in ( 12 ) with (/(X^,) replaced by 
exp(— g 2 De) and exp {—r] 2 D G ) respectively. The threshold 
finding algorithm in (19) is altered accordingly. 

IV. Mechanizing the IS algorithm 
A. A methodology 

The various issues of the preceding sections are summarized 
in an IS simulation methodology that outlines the principal 
steps required for implementation of the adaptive algorithms. 
It is a cautious methodology which has been used in previous 
works on applications of IS. 

Of interest is the extreme but realistic (and often encoun¬ 
tered) situation where we have no knowledge whatsoever of 
the FAP a of a detection algorithm for a given value of thresh¬ 
old ?/. In such a case, referring to the objective function of 
Figure 1, the basic idea of fast adaptive simulation is to travel 
down the curve from its left side. An initial (a 0 ,r] 0 ) pair is 
needed. It is easily obtained using conventional MC simulation 
for a high value of a 0 , say 0.1 (and a correspondingly low 
value of r ] 0 ), using K = 100/a o = 1000 trials according to 
the well known MC thumb rule. This can be accomplished 
quickly with a few experiments. The (a 0 ,rj 0 ) pair (whose 
accuracy need not be very high) provides the starting values 
for the IS procedure. It begins by forming the estimate I g (9) at 
the threshold ?/ G (using 1000 trials) and locating its minimum 
(possibly graphically) as a function of the biasing parameter 9. 
It is advisable to search for the optimum bias, denoted as 9 0 , 
starting from 9 = 1 to avoid locating an overbiasing minimizer. 
This leads, via (23), to an estimate 1 of the maximum IS 
gain available (with a possible correction for a) and results 
in the quadruplet (a 0 ,r] 0 ,9 0 ,T(9 0 )). At this stage the IS 
adaptations are implemented. The threshold-finding algorithm 
of (19) is implemented with initial value r] 0 for a new pre¬ 
specified a 0 slightly smaller than the initial a 0 and using 
a conservative number of K = lOO/(a o r(0 o )) for the IS 
trials. It is conservative because we know that simulation gain 
increases with decreasing rare-event probability and hence this 
number is guaranteed to provide better than the rule-of-thumb 
accuracy at the new a 0 . Note that most of the IS gain can be 
leveraged if the FAP decrements are kept small. The biasing 
algorithm of (17) is implemented simultaneously using 9 0 as 
initial value. At the end of the adaptations the resulting gain 
is estimated for updating K and a new quadruplet is obtained. 
This is continued until we have a complete characterization of 
false-alarm behaviour of the detector down to the desired a 0 . 
The procedure is summarized below. 

( \ p 

1) Implementation: Define the set {cto } 1 such that 0.1 = 

a'o' > > > • • • > a'J 3 ' 1 = 10 ~ 7 , for example. 

Pre-processing: 

Step 1. p = 1. Use 1000 conventional MC trials to 

obtain (a^, pair. 

Step 2. Find 9^ = argmin I g (9) and calculate 

6 

?(p\ 9 { 0 p) ) using (23). 

Adaptive IS: 

Step 3. p = p + 1 
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Step 4. Let m = 1. Set 9^ = 9o P 1 \ = po’ ^ 

and A'= 100/(ai p) f(P- 1 )(6'i P_1) )). 

Step 5. Generate K secondary vectors and compute 
B and I) using (8) and (13). 

Step 6. Implement (12), (18) using biased secondary 
vectors with parameter 9 $. 

(p) (p) 

Step 7. Compute 9£ +l and tj^' +1 using (17) and 

(p) 

(19); if both algorithms have converged, set 9 0 = 

<&. vP = Cli- r ( p\ei p) ) = f(p)(0^ 1 ) using 

(23), and go to step 3 or stop if p = P\ otherwise 
go to step 5 with m = m + 1. 

At the end of simulations we have the P quadruplets 

{ai P \r,i P \9 { 0 p \f(p)(9i p) )} P 1 . 


B. Complexity reduction for IS adaptations 

The methodology is simple enough but certain observations 
can be made which lead to an artifice that further reduces 
computational effort to a considerable extent. In rare-event 
simulation of highly reliable systems that involve complex 
signal processing operations, two factors contribute to simula¬ 
tion time. The first is the rare event itself that is under study; 
this is handled by suitable IS biasing techniques. The second 
factor is the computational intensity of the signal processing 
of the system. In STAP detectors the chief processing burden 
is from inversion of large matrices. Several millions of MC 
trials with as many matrix inversions are needed to estimate 
low FAPs. Using an IS scheme the number of trials can 
be reduced appreciably, for a single rare-event probability 
estimation. However, there remains the important matter of 
executing the IS adaptations of (17) and (19), using this 
small number of IS trials for each iteration. These can be 
computationally burdensome if the recursions are to be carried 
out for several a 0 values. This is where the artifice comes in. 
It turns out, in contrast to conventional MC simulation, that 
adaptations can be performed by reusing the unbiased random 
quantities generated for any one set of trials. Consequently it is 
not necessary to run truly randomized adaptive IS algorithms. 
The result will be a large savings in computational effort. 

This idea is illustrated by rewriting the FAP estimator of 
(12) in the form 



i=l 


- / (25) 


obtained by substituting (14) into (12), replacing D and B 
respectively by their unbiased versions denoted as D v and 
B v , and observing that whereas the estimator in (12) uses 
samples drawn from biased distributions /*, the mathemati¬ 
cally equivalent estimator of (25) operates with only unbiased 
quantities, obtained via /. Biasing is handled by the (explicit) 
presence of 9. The number of trials remains unchanged, at 
K. The estimators required in (17) and (19) can be rewritten 
in a similar manner. It is a straightforward exercise and is 
omitted. Assuming then that the value of K is large enough 
for estimation of the target a 0 (and this is guaranteed by 
the methodology described above), it follows that the set 
of K instances of D v and B v employed in (25) can be 


repeatedly used in the adaptations. With the complexity of 
frequent regeneration thus removed, these adaptations (such 
as minimization of I g (9) in (18)) can be extremely fast, 
often requiring only a few iterations. It is assumed that a 
sufficiently large pool of unbiased variables is pre-generated 
to accommodate K trials. 

There are two related issues that are briefly discussed here 
in qualitative terms. In IS simulations carried out here (and 
elsewhere), it is tacitly assumed that the number of trials I\ is 
also large enough to accurately estimate I g (9). This quantity is 
admittedly not a rare-event probability and its estimator I g (9), 
given in (18), is clearly unbiased. Nevertheless, its accuracy 
should be studied in terms of its error variance, irrespective 
of whether we employ the above described reuse strategy or 
not. In this regard it must be pointed out that the biasing 
densities used for IS simulations are being chosen to minimize 
I g (9) (or its estimate I g (9)), and not to estimate it accurately. 
Indeed, attempting to estimate I g (9) accurately would lead us 
to consider further biasing of the original biasing distribution. 
As probably known to IS researchers, the more powerful a 
biasing distribution is (in the sense of being “close” to the 
unrealizable optimal biasing density that estimates the rare- 
event probability exactly), the more “constant” will be the 
i.i.d. random variables yi defined after (24). Interestingly this 
implies that even for small K, the estimated variance var a g 
can become small for such biasing distributions. The latter 
observations accentuate the need to search for good simulation 
distributions. Unfortunately, it is beyond the remit of this paper 
to investigate the variance and rate of convergence of I g (9) and 
the dependencies of these on biasing distributions and number 
of trials K. Some remarks on this issue can be found in [16]. 
We have confirmed the legitimacy of the reuse strategy in 
several cases in two ways. Using an overkill, that is, employing 
a large number of IS trials for varying 9 , it has been ascertained 
that the correct value of the optimum bias does not differ 
noticeably from the 9 0 determined as above. Secondly, truly 
randomized optimization algorithms have been used, with the 
same conclusion. 

The second issue concerns the ease with which the form of 
the estimator in (25) is obtained. The reason lies in the biasing 
technique used. Biasing in IS can be classed, for the purpose 
of this discussion, into two types. One kind derives biased 
quantities by direct parameterized transformations of unbiased 
variables, such as the scaling employed here. The second is 
just the complementary class, where such transformations do 
not exist. Examples of the latter can include new distributions 
for conducting simulations and some instances of distributions 
based on large deviations theory arguments. Within the first 
class, only the subset of transformations that are many-to-one 
will allow placing the IS estimator in a form such as that of 
(25). 


V. Simulation results 

Implementation results of the fast simulation procedures 
using the //-method are described here. A typical example 
of evolution of the threshold-finding algorithm is shown in 
Fig. 2 for the square-law AMF detector with L = 704 training 
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vectors and space-time product N = 352. Clearly, about 20 
recursions of (19) appear sufficient for convergence. Note also 
from this figure that the “converged” value of the threshold 
for an a 0 acts as the initial value for the next lower a 0 , as 
described in Section IV-A. Although not shown here, similar 
estimation results have been obtained for the E-AMF and 
GM-STAP detectors. The results for the AMF detector almost 
“coincide” with numerical computations of (2). As analytical 
expressions or approximations for the FAPs are not available, 
such a comparison is not possible for the variants. Threshold 
multiplier values for the 3 detectors are summarized in Fig. 3 
for L = 128 and N = 64. For the E-AMF and GM-STAP 
detectors the interpolations 

rj E = 0.0095 x 3 - 0.1713 a ; 2 + 1.8869 x + 1.6873 

rj G = 0.0111 x 3 - 0.1985 a ; 2 + 2.2351 x + 1.9899 

for L = 128 and N = 64, and 

rj E = 0.01a : 3 - 0.1865 a ; 2 + 1.8822 x + 1.7255 

i] G = 0.0118 a ; 3 - 0.2198 a ; 2 + 2.2292 x + 2.0408 

for L = 704 and N = 352, obtained from inverse IS results 
where x = — log 10 a, can be used for estimating multipliers 
for FAPs in the range 10 _1 to 10 -7 . 

Performance graphs of the IS algorithms are in Figures 4 
and 5, which show the estimated IS gain T obtained over 
MC simulation and the number of trials required respectively. 
The number of IS trials, although appreciably reduced in 
comparison with MC, increases rapidly with decrease in false 
alarm probability. From these two figures it can be observed 
that, for a FAP of 10 -6 , the square-law AMF with L = 128 
and IV = 64 requires 1500 IS trials whereas only 100 trials are 
needed for L = 704 and N = 352. This can be attributed to 
the fact that D and B, the only random quantities that appear 
in ( 12 ), have more concentrated density functions in the latter 
case, thus causing the IS biasing to be more effective during 
simulation. 

Detection probabilities for all three detectors have been 
estimated and compared at a FAP of 10 _B . For homogeneous 
clutter backgrounds and non-fluctuating (Swerling 0) and 
fluctuating (Swerling 1) targets, the results are in Fig. 6 . The 
SNR loss of the two variants compared to the AMF detector 
is very slight, the maximum being 0.3 db for a Swerling 0 
target. Performance degradations for training data contami¬ 
nated by nonhomogeneities consisting of interfering targets 
that are assumed to have the same Doppler-angle properties 
and characteristics as the primary target are shown in Fig. 7 
for two Swerling 0 and Swerling 1 targets. Each interferer is 
assumed to have the same steering vector and power as the 
primary target. Similar results (not shown here) have been 
obtained for different numbers of interferers. It is evident 
that the GM-STAP detector is most robust in the presence of 
interferers, enjoying (in some cases) several decibels of power 
advantage over the square-law AMF detector. Consequently, 
its FAP performance in these situations (though not evaluated 
here) should also be relatively robust. 


A. Comments 

Despite the accuracy of our results, biasing by scaling has 
not done as well for the variants as for the AMF detector. 
Furthermore, even for the latter, IS performance should be im¬ 
proved by a better choice of biasing scheme in order to achieve 
the goal of having the required number of trials to be constant, 
or at least increasing very slowly, with decrease of FAP. 
This is certainly a matter for further investigation that may 
pay dividends while examining other detector configurations. 
Scaling was used in this work in the interest of expediency. 
Some details regarding the behaviour of the scaling parameter 
6 are available in [15], Briefly, it is close to unity and has a 
small spread over the range of FAPs considered. This is due 
to the shape of the density functions of B and D. 

VI. Conclusion 

A humble inroad into the use of adaptive IS algorithms 
to characterize a STAP detector has been made. The AMF 
detector was used for validation and results have been good. 
The chief reasons for this are that we were able to invoke the 
(/-method and inverse IS, use a biasing strategy that could 
be easily optimized adaptively, and find a way around the 
difficult task of inverting large matrices several times during 
simulations. As a small demonstration of the potential of IS, 
two AMF detector variants that are known to be relatively 
obdurate to mathematical analysis have been suggested and 
characterized. Sufficient numerical results have been provided 
here in order to motivate interested researchers to confirm 
our findings and possibly carry out their own investigations 
into the subject. There are some technical issues regarding IS 
estimators that remain to be settled; however, the development 
of better biasing methods is an important matter. In any case, 
our hope is that application of these fast simulation techniques 
to more advanced STAP configurations will also meet with 
success. But this remains to be seen as we are certainly not in 
position to predict what subtleties (and difficulties) these other 
detection algorithms can throw up. It is clear that IS is still 
in its infancy, especially insofar as its use for characterizing 
modern detection algorithms is concerned. 

Appendix I 
CFAR PROPERTY 

An invariance property is established that can be used to 
show that certain STAP detection algorithms have FAPs that 
do not depend on the data covariance R. The proposition 
given here follows quite simply from arguments contained 
in the exposition of the generalized likelihood ratio STAP 
detector test (Kelly’s GLRT) given in [1], They are outlined 
here for convenience of the reader 2 . Assume, as before, that the 
primary and training data vectors have the same covariance. 
Consider the variables 

G = si'RT^x and G(l) = s t R" 1 x(/) (26) 

2 It would be helpful for the reader to refer to [1] (see also [2]). We have used 
several results from this now classic paper and have attempted to maintain 
the same notation. 
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for l = 1Using the transformations u = R 1 / 2 s, 
y = R~ 1 / 2 x, and y(Z) = R _ 1 ( 2 x(Z), leads to 

G = u^R' 1 y and G(Z) = u t R" 1 y(0 (27) 


where R = R'7 2 RR _1 ( 2 . The whitened vectors Y and 
Y(Z) are both distributed CAf^O, I). It turns out that R 
has the complex Wishart distribution CW(L,N;^ I), [13]. 
Further, a unitary transformation U can be found which 
rotates the new signal vector u into an elementary vector 
e as de = U^u, such that e = [1,0 ,..., 0] 7 and where 
d 2 = u u = s'R As. The first column of U is the new signal 
vector u/d. The remaining columns comprise an orthonormal 
basis determined, for example, by a Gram-Schmidt procedure. 
Let z = lAy and z(Z) = U^y(0- Applying these to (27) 
yields the variables 

G=ye J< S^ 1 z and G(Z) = y e^S~ 1 z(l) (28) 

1 j Lj 


where S = L U'RU. While Z and Z(Z) are distributed 
as CAf/yfO. I) and are independent, S has the distribution 
CW(L,N\1). The vectors z and z(l) are decomposed as 
z = [ z A z b ] t and z(Z) = [z 4 (Z) z B (l)] J where the A 
components are scalar and B components (N — 1)-vector. 
Correspondingly, S is decomposed as 


s = E z W z W f = 


i=i 


^AA $ab 

c c 


(29) 


so that S AB = Ef =1 Z A (0 Z B (0 f > S BB = Zti Z B (0 z b (0 f 
and so on. Also 


e's- 1 = [r AA v AB ] 

where 

^ AA ~ (^AA — ^AB^BB^Ba) 
'PAB ~ ^AA^AB^BB 


(30) 


(31) 


and 

E B {\Y(l)\ 2 }=l-z B (y(Y / z B (i)z B (i)A \ B (l) 

'»=.l ' 

for l = 1,... ,L, with E g denoting conditional expectation. 
Further, the conditional covariance of the variables Y(l) is 
given by 


E B {Y(k)Y(nY} 


1=1 


z B (fc) 


for k Y n. 

Hence the set of conditionally jointly Gaussian zero mean 

random variables Y and {Y(Z)} have individual variances 

and covariances that are functions of the random vectors z D 
L B 
and {z B (/)} 1 . The latter are all jointly independent, each 

being distributed as CA/"jv-i(0,1). The probability of any 
event defined on the random variables Y and {Y (Z)} in (33) 
can thus be determined by performing an averaging operation 
over the distributions of z B and {z B (Z)} and this probability 
will be independent of the data covariance R. This statement 
is also true for the random variables G and {G(Z)} in (32) 
with the caveat that any constant scaling of these variables 
should leave the event unchanged. The preceding arguments 
therefore constitute proof of the following 
Proposition: Any STAP detection algorithm that uses only 
the random variables G and {G(Z)| defined in (26) for its 
description such that the algorithm itself is unchanged by 
arbitrary but equal scaling of all these variables, has a FAP 
which is independent of the target-free data covariance R. 

It follows immediately from this that both the E-AMF and 
GM-STAP detectors have FAPs that are independent of the 
covariance matrix R. 


Appendix II 

Adaptive algorithms for 2-d biasing 


which follow from the Frobenius relations for inversion of 
block matrices, [1], Substituting (30) and (31) in (28) gives 

G=d L V A A y and G{l)= d L V AA y{l) (32) 

where 

L 

y = Z A-^2 Z A( l > B ( l ) fS BB Z B 
1=1 

(33) 

L 

y(l) = z A {l)-Y J z A^ z B^ S B l B z B {l) 
i-1 

Conditioned on the vectors z B and z B (l), it follows that the 
random variables Y and Y(l) in (33) are (in the absence of 
target) zero mean Gaussian. With a little more algebra it can 
be shown, as in [ 1 ], that they are uncorrelated with variances 


From the first line of (10), the /-function is estimated as 
1 K 

= < 34 > 

i=1 

and its minimization can be carried out using the 2 - 
dimensional adaptive algorithm 

*Wi = v m - <5J“ 1 V/(iz m ), m= 1,2,... (35) 

Here, S is a step-size parameter used to control convergence, 
and to is the index of iteration. This is a stochastic Newton 
recursion. It achieves minimization of I by estimating a 
solution of 

vj(i/) = (/„ T e ) T = o 

where I a = dl(v)/da and Ig = dl(v)/d0. The estimate of 
the Jacobian J (which is the Hessian matrix of I) is given by 


E B m 2 } 


1 + z 


L 

B \E Z B 

1=1 


(i) z B oy 


j = 
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where I xy = dl x /dy. Starting with the second line of (10), 
the various /-functions defined above can be obtained by using 
the notational equations 

4 = E{l(A)dW/dx} 

= E+{l{A)WW x } 

and similarly 

Ixx = E*{1(A)WW XX } 

Ixy = E*{l(A)WW X y} 

with the corresponding weighting function derivatives easily 
calculated, using (7), as W x = dW/dx, W xy = d 2 W/dxdy, 
and so on. These /-functions can be estimated as in (34). 

Appendix III 

Gain of the ^-method estimator 
The proof given here is in the context of STAP detectors 
and is a generalization of the one in [10] for conventional 
CFAR detection algorithms. A more complete generalization 
for estimators in an arbitrary rare-event setting is available in 
[11]. 

When no IS is used (W = 1), then I g = E{g 2 (X.L)} < 
E{g(X.L)} = a. To show that I g < I = I(u) (the biasing 
vector is omitted for convenience) with IS, some additional 
notation is useful. Let U = | s' R~ 1 X | 2 and let V denote the 
right hand side of any of the three STAP detectors, excluding 
the multiplier. We bear in mind that U is a function of X and 
X/,, and V a function of X/ . Then, in the absence of target, 
a = P(U > ijV) and (10) can be written as 

I = E{l(U>r}V)W(X,X L )} 

with W defined in the first line of (7). For the (/-method 
estimator, we can write 

g(x L ) = P(U > yV | X L = x L ) 

= E{l(U>r/V) IX L = x L } 

= J l(u>r]v)f{x\x L )dx 

= / l(u > gv) W c (x,x L ) /*(x|x L ) dx 

where W c (x,x L ) = /(x|x L )//*(x|x L ). Then 

5 2 ( x l) < f 1(« > ilv) W c 2 (x,x L )/*(x|x L ) dx 

= J l(w > T)v) W c (x,x L ) /(x|x L ) dx 

by Jensen’s inequality. Therefore I g , defined in (16), is given 
by 

!g = J 3 2 (xl) W(x l ) f(x L ) dx L 

<-Ii 
-II 
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1 (u > rjv) W c (x,x L ) W(x L ) f(x |x L ) 

• f(x L ) dxdx L 

1 (u > r)v) W (x, Xi) f(x,x L ) dxdx L 



IEEE TRANSACTIONS ON AES. TO APPEAR IN JANUARY 2007 10 

List of Figures 

1 Example objective function (for the E-AMF detector). 11 

2 Threshold optimization for AMF detector (L = 704, N = 352) using inverse IS. 11 

3 Multipliers for L = 128 and N = 64. The starred dots indicate validations for the E-AMF detector using 

conventional MC simulations. 12 

4 IS gains. The set of 3 lower graphs correspond to L = 128 and TV = 64. 12 

5 Simulation lengths for L = 128 and N = 64. 13 

6 Detection performances in homogeneous clutter. L = 128, N = 64, and a 0 = 10 6 . 13 

7 Detection with 2 interfering targets in training data. L = 128, N = 64, and a Q = 10 6 . 14 










IEEE TRANSACTIONS ON AES, TO APPEAR IN JANUARY 2007 


x 10 12 



Fig. 1. Example objective function (for the E-AMF detector). 


AMF: L = 704, N = 352 



Fig. 2. Threshold optimization for AMF detector (L = 704, AT = 352) using inverse IS. 
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Fig. 3. Multipliers for L = 128 and N = 64. The starred dots indicate validations for the E-AMF detector using conventional MC simulations. 



Fig. 4. IS gains. The set of 3 lower graphs correspond to L = 128 and N = 64. 
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Fig. 5. Simulation lengths for L = 128 and N = 64. 



Fig. 6. Detection performances in homogeneous clutter. L = 128, N = 64, and a Q = 10 6 . 
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Fig. 7. Detection with 2 interfering targets in training data. L = 128, N = 64, and a Q = 10 6 . 




















































